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: Oct 02 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){ not mapping yet */
 143         if(pch==ech){
 144
 145           /* we need a 1/4 envelope window overlap front and back */
 146           double *pcm=v->pcm[pch]+dcurr*step-step/2;
 147           _ve_deltas(mult,pcm,dtotal-dcurr,window,winlen);
 148
 149         }
 150       }
 151     }
 152     v->envelope_current=dtotal;
 153     frame++;
 154   }
 155 }
 156
 157 /* This readies the multiplier vector for use/coding.  Clamp/adjust
 158    the multipliers to the allowed range and eliminate unneeded
 159    coefficients */
 160
 161 void _ve_envelope_sparsify(vorbis_block *vb){
 162   int ch;
 163   for(ch=0;ch<vb->vd->vi.envelopech;ch++){
 164     int flag=0;
 165     double *mult=vb->mult[ch];
 166     int n=vb->multend;
 167     double first=mult[0];
 168     double last=first;
 169     double clamp;
 170     int i;
 171
 172     /* are we going to multiply anything? */
 173
 174     for(i=1;i<n;i++){
 175       if(mult[i]>=last*vb->vd->vi.preecho_thresh){
 176         flag=1;
 177         break;
 178       }
 179       if(i<n-1 && mult[i+1]>=last*vb->vd->vi.preecho_thresh){
 180         flag=1;
 181         break;
 182       }
 183       last=mult[i];
 184     }
 185
 186     if(flag){
 187       /* we need to adjust, so we might as well go nuts */
 188
 189       int begin=i;
 190       clamp=last?last:1;
 191
 192       for(i=0;i<begin;i++)mult[i]=0;
 193
 194       last=1;
 195       for(;i<n;i++){
 196         if(mult[i]/last>clamp*vb->vd->vi.preecho_thresh){
 197           last=mult[i]/vb->vd->vi.preecho_clamp;
 198
 199           mult[i]=floor(log(mult[i]/clamp/vb->vd->vi.preecho_clamp)/log(2))-1;
 200           if(mult[i]>15)mult[i]=15;
 201         }else{
 202           mult[i]=0;
 203         }
 204       }
 205     }else
 206       memset(mult,0,sizeof(double)*n);
 207   }
 208 }
 209
 210 void _ve_envelope_apply(vorbis_block *vb,int multp){
 211   vorbis_info *vi=&vb->vd->vi;
 212   double env[vb->multend*vi->envelopesa];
 213   envelope_lookup *look=&vb->vd->ve;
 214   int i,j,k;
 215
 216   for(i=0;i<vi->envelopech;i++){
 217     double *mult=vb->mult[i];
 218     double last=1.;
 219
 220     /* fill in the multiplier placeholders */
 221
 222     for(j=0;j<vb->multend;j++){
 223       if(mult[j]){
 224         last=mult[j];
 225       }else
 226         mult[j]=last;
 227     }
 228
 229     /* compute the envelope curve */
 230     if(_ve_envelope_generate(mult,env,look->window,vb->multend,
 231                              vi->envelopesa)){
 232
 233       /* apply the envelope curve */
 234       for(j=0;j<vi->channels;j++){
 235
 236         /* check to see if the generated envelope applies to this channel */
 237         /*if(vi->envelopemap[j]==i){ not mapping yet */
 238         if(j==i){
 239
 240           if(multp)
 241             for(k=0;k<vb->multend*vi->envelopesa;k++)
 242               vb->pcm[j][k]*=env[k];
 243           else
 244             for(k=0;k<vb->multend*vi->envelopesa;k++)
 245               vb->pcm[j][k]/=env[k];
 246         }
 247       }
 248     }
 249   }
 250 }
 251