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 07 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 void _ve_envelope_clear(envelope_lookup *e){
  47   if(e->window)free(e->window);
  48   memset(e,0,sizeof(envelope_lookup));
  49 }
  50
  51 /* initial and final blocks are special cases. Eg:
  52    ______
  53          `--_
  54    |_______|_`-.___|_______|_______|
  55
  56               ___
  57           _--'   `--_
  58    |___.-'_|_______|_`-.___|_______|
  59
  60                       ___
  61                   _--'   `--_
  62    |_______|___.-'_|_______|_`-.___|
  63
  64                                _____
  65                            _--'
  66    |_______|_______|____.-'|_______|
  67
  68    as we go block by block, we watch the collective metrics span. If we
  69    span the threshhold (assuming the threshhold is active), we use an
  70    abbreviated vector */
  71
  72 static int _ve_envelope_generate(double *mult,double *env,double *look,
  73                                  int n,int step){
  74   int i,j,p;
  75   double m,mo;
  76
  77   for(j=0;j<n;j++)if(mult[j]!=1)break;
  78   if(j==n)return(0);
  79
  80   n*=step;
  81   /* first multiplier special case */
  82   m=ldexp(1,mult[0]-1);
  83   for(p=0;p<step/2;p++)env[p]=m;
  84
  85   /* mid multipliers normal case */
  86   for(j=1;p<n-step/2;j++){
  87     mo=m;
  88     m=ldexp(1,mult[j]-1);
  89     if(mo==m)
  90       for(i=0;i<step;i++,p++)env[p]=m;
  91     else
  92       for(i=0;i<step;i++,p++)env[p]=m*look[i]+mo*look[i+step];
  93   }
  94
  95   /* last multiplier special case */
  96   for(;p<n;p++)env[p]=m;
  97   return(1);
  98 }
  99
 100 /* right now, we do things simple and dirty (read: our current preecho
 101    is a joke).  Should this prove inadequate, then we'll think of
 102    something different.  The details of the encoding format do not
 103    depend on the exact behavior, only the format of the bits that come
 104    out.
 105
 106    Mark Taylor probably has much witter ways of doing this...  Let's
 107    see if simple delta analysis gives us acceptible results for now.  */
 108
 109 static void _ve_deltas(double *deltas,double *pcm,int n,double *win,
 110                        int winsize){
 111   int i,j,p=winsize/2;
 112   for(j=0;j<n;j++){
 113     p-=winsize/2;
 114     for(i=0;i<winsize-1;i++,p++){
 115       double temp=fabs(win[i]*pcm[p]-win[i+1]*pcm[p+1]);
 116       if(deltas[j]<temp)deltas[j]=temp;
 117     }
 118     p++;
 119   }
 120 }
 121
 122 void _ve_envelope_multipliers(vorbis_dsp_state *v){
 123   int step=v->samples_per_envelope_step;
 124   static int frame=0;
 125
 126   /* we need a 1-1/4 envelope window overlap begin and 1/4 end */
 127   int dtotal=(v->pcm_current-step/2)/v->samples_per_envelope_step;
 128   int dcurr=v->envelope_current;
 129   double *window=v->ve.window;
 130   int winlen=v->ve.winlen;
 131   int pch,ech;
 132   vorbis_info *vi=&v->vi;
 133
 134   if(dtotal>dcurr){
 135     for(ech=0;ech<vi->envelopech;ech++){
 136       double *mult=v->multipliers[ech]+dcurr;
 137       memset(mult,0,sizeof(double)*(dtotal-dcurr));
 138
 139       for(pch=0;pch<vi->channels;pch++){
 140
 141         /* does this channel contribute to the envelope analysis */
 142         if(vi->envelopemap[pch]==ech){
 143
 144           /* we need a 1/4 envelope window overlap front and back */
 145           double *pcm=v->pcm[pch]+dcurr*step-step/2;
 146           _ve_deltas(mult,pcm,dtotal-dcurr,window,winlen);
 147
 148         }
 149       }
 150     }
 151     v->envelope_current=dtotal;
 152     frame++;
 153   }
 154 }
 155
 156 /* This readies the multiplier vector for use/coding.  Clamp/adjust
 157    the multipliers to the allowed range and eliminate unneeded
 158    coefficients */
 159
 160 void _ve_envelope_sparsify(vorbis_block *vb){
 161   int ch;
 162   for(ch=0;ch<vb->vd->vi.envelopech;ch++){
 163     int flag=0;
 164     double *mult=vb->mult[ch];
 165     int n=vb->multend;
 166     double first=mult[0];
 167     double last=first;
 168     double clamp;
 169     int i;
 170
 171     /* are we going to multiply anything? */
 172
 173     for(i=1;i<n;i++){
 174       if(mult[i]>=last*vb->vd->vi.preecho_thresh){
 175         flag=1;
 176         break;
 177       }
 178       if(i<n-1 && mult[i+1]>=last*vb->vd->vi.preecho_thresh){
 179         flag=1;
 180         break;
 181       }
 182       last=mult[i];
 183     }
 184
 185     if(flag){
 186       /* we need to adjust, so we might as well go nuts */
 187
 188       int begin=i;
 189       clamp=last?last:1;
 190
 191       for(i=0;i<begin;i++)mult[i]=0;
 192
 193       last=1;
 194       for(;i<n;i++){
 195         if(mult[i]/last>clamp*vb->vd->vi.preecho_thresh){
 196           last=mult[i]/vb->vd->vi.preecho_clamp;
 197
 198           mult[i]=floor(log(mult[i]/clamp/vb->vd->vi.preecho_clamp)/log(2))-1;
 199           if(mult[i]>15)mult[i]=15;
 200         }else{
 201           mult[i]=0;
 202         }
 203       }
 204     }else
 205       memset(mult,0,sizeof(double)*n);
 206   }
 207 }
 208
 209 void _ve_envelope_apply(vorbis_block *vb,int multp){
 210   vorbis_info *vi=&vb->vd->vi;
 211   double env[vb->multend*vi->envelopesa];
 212   envelope_lookup *look=&vb->vd->ve;
 213   int i,j,k;
 214
 215   for(i=0;i<vi->envelopech;i++){
 216     double *mult=vb->mult[i];
 217     double last=1.;
 218
 219     /* fill in the multiplier placeholders */
 220
 221     for(j=0;j<vb->multend;j++){
 222       if(mult[j]){
 223         last=mult[j];
 224       }else
 225         mult[j]=last;
 226     }
 227
 228     /* compute the envelope curve */
 229     if(_ve_envelope_generate(mult,env,look->window,vb->multend,
 230                              vi->envelopesa)){
 231
 232       /* apply the envelope curve */
 233       for(j=0;j<vi->channels;j++){
 234
 235         /* check to see if the generated envelope applies to this channel */
 236         if(vi->envelopemap[j]==i){
 237
 238           if(multp)
 239             for(k=0;k<vb->multend*vi->envelopesa;k++)
 240               vb->pcm[j][k]*=env[k];
 241           else
 242             for(k=0;k<vb->multend*vi->envelopesa;k++)
 243               vb->pcm[j][k]/=env[k];
 244         }
 245       }
 246     }
 247   }
 248 }
 249