postfish/deverb.c

   1 /*
   2  *
   3  *  postfish
   4  *
   5  *      Copyright (C) 2002-2005 Monty and Xiph.Org
   6  *
   7  *  Postfish is free software; you can redistribute it and/or modify
   8  *  it under the terms of the GNU General Public License as published by
   9  *  the Free Software Foundation; either version 2, or (at your option)
  10  *  any later version.
  11  *
  12  *  Postfish is distributed in the hope that it will be useful,
  13  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
  14  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
  15  *  GNU General Public License for more details.
  16  *
  17  *  You should have received a copy of the GNU General Public License
  18  *  along with Postfish; see the file COPYING.  If not, write to the
  19  *  Free Software Foundation, 675 Mass Ave, Cambridge, MA 02139, USA.
  20  *
  21  *
  22  */
  23
  24 #include "postfish.h"
  25 #include "feedback.h"
  26 #include <fftw3.h>
  27 #include "subband.h"
  28 #include "bessel.h"
  29 #include "deverb.h"
  30
  31 /* (since this one is kinda unique) The Deverberation Filter....
  32
  33    Reverberation in a measurably live environment displays
  34    log amplitude decay with time (linear decay when plotted on a dB
  35    scale).
  36
  37    In its simplest form, the deverber follows actual RMS amplitude
  38    attacks but chooses a slower-than-actual decay, then expands
  39    according to the dB distance between the slow and actual decay.
  40
  41    Thus, the deverber can be used to 'dry out' a very 'wet'
  42    reverberative track. */
  43
  44 extern int input_size;
  45 extern int input_rate;
  46 extern int input_ch;
  47
  48 typedef struct {
  49   subband_state ss;
  50
  51   iir_filter smooth;
  52   iir_filter smoothlimit;
  53   iir_filter release;
  54
  55   iir_state *iirS[deverb_freqs];
  56   iir_state *iirR[deverb_freqs];
  57
  58   float prevratio[deverb_freqs];
  59
  60   int inactive_delay[deverb_freqs];
  61
  62 } deverb_state;
  63
  64 deverb_settings deverb_channel_set;
  65 static deverb_state channel_state;
  66 static subband_window sw;
  67
  68 static void deverb_reset_band(int freq){
  69   int j;
  70   for(j=0;j<input_ch;j++)
  71     memset(&channel_state.iirS[freq][j],0,sizeof(iir_state));
  72   for(j=0;j<input_ch;j++)
  73     memset(&channel_state.iirR[freq][j],0,sizeof(iir_state));
  74   channel_state.inactive_delay[freq]=2;
  75 }
  76
  77 static void deverb_reset_one(int freq,int ch){
  78   memset(&channel_state.iirS[freq][ch],0,sizeof(iir_state));
  79   memset(&channel_state.iirR[freq][ch],0,sizeof(iir_state));
  80 }
  81
  82 void deverb_reset(){
  83   int i;
  84
  85   subband_reset(&channel_state.ss);
  86
  87   for(i=0;i<deverb_freqs;i++)
  88       deverb_reset_band(i);
  89 }
  90
  91 static void filter_set(subband_state *ss,
  92                        float msec,
  93                        iir_filter *filter,
  94                        int attackp,
  95                        int order){
  96   float alpha;
  97   float corner_freq= 500./msec;
  98
  99   /* make sure the chosen frequency doesn't require a lookahead
 100      greater than what's available */
 101   if(impulse_freq2(input_size*2-ss->qblocksize*3)*1.01>corner_freq &&
 102      attackp)
 103     corner_freq=impulse_freq2(input_size*2-ss->qblocksize*3);
 104
 105   alpha=corner_freq/input_rate;
 106   filter->g=mkbessel(alpha,order,filter->c);
 107   filter->alpha=alpha;
 108   filter->Hz=alpha*input_rate;
 109   filter->ms=msec;
 110 }
 111
 112 int deverb_load(void){
 113   int i;
 114   int qblocksize=input_size/16;
 115   memset(&channel_state,0,sizeof(channel_state));
 116
 117   deverb_channel_set.active=calloc(input_ch,sizeof(*deverb_channel_set.active));
 118
 119   subband_load(&channel_state.ss,deverb_freqs,qblocksize,input_ch);
 120   subband_load_freqs(&channel_state.ss,&sw,deverb_freq_list,deverb_freqs);
 121
 122   for(i=0;i<deverb_freqs;i++){
 123     channel_state.iirS[i]=calloc(input_ch,sizeof(iir_state));
 124     channel_state.iirR[i]=calloc(input_ch,sizeof(iir_state));
 125     channel_state.inactive_delay[i]=2;
 126   }
 127   return 0;
 128 }
 129
 130 static void deverb_work_helper(void *vs, deverb_settings *sset){
 131   deverb_state *sss=(deverb_state *)vs;
 132   subband_state *ss=&sss->ss;
 133   int i,j,k,l;
 134   float smoothms=sset->smooth*.1;
 135   float releasems=sset->release*.1;
 136   iir_filter *smooth=&sss->smooth;
 137   iir_filter *smoothlimit=&sss->smoothlimit;
 138   iir_filter *release=&sss->release;
 139   int ahead;
 140
 141   if(smoothms!=smooth->ms){
 142     filter_set(ss,smoothms,smooth,1,2);
 143     filter_set(ss,smoothms,smoothlimit,0,1);
 144   }
 145   if(releasems!=release->ms)filter_set(ss,releasems,release,0,1);
 146
 147   ahead=impulse_ahead2(smooth->alpha);
 148
 149   for(i=0;i<deverb_freqs;i++){
 150     int firstlink=0;
 151     float fast[input_size];
 152     float slow[input_size];
 153     int ratio=sset->ratio[i];
 154     float multiplier = (sss->inactive_delay[i]>0 ? 0. : 1.-1000./ratio);
 155
 156     if(ratio==1000 && sss->prevratio[i]==0.){
 157
 158       /* although the effect is active, this band is set to unity.
 159          Save CPU: Do not run the filters.  This is not tied into the
 160          subbanding 'active' indicator as that's by channel, not by
 161          band. */
 162
 163       deverb_reset_band(i);
 164
 165     }else{
 166       sss->inactive_delay[i]--;
 167       if(sss->inactive_delay[i]<0)sss->inactive_delay[i]=0;
 168
 169       for(j=0;j<input_ch;j++){
 170         int active=(ss->effect_active1[j] ||
 171                     ss->effect_active0[j] ||
 172                     ss->effect_activeC[j]);
 173
 174         if(!active){
 175           deverb_reset_one(i,j);
 176         }else{
 177
 178           /* run the filters */
 179
 180           if(sset->linkp){
 181             if(!firstlink){
 182               firstlink++;
 183               memset(fast,0,sizeof(fast));
 184               float scale=1./input_ch;
 185               for(l=0;l<input_ch;l++){
 186                 float *x=sss->ss.lap[i][l]+ahead;
 187                 for(k=0;k<input_size;k++)
 188                   fast[k]+=x[k]*x[k];
 189               }
 190               for(k=0;k<input_size;k++)
 191                 fast[k]*=scale;
 192
 193             }
 194
 195           }else{
 196             float *x=sss->ss.lap[i][j]+ahead;
 197             for(k=0;k<input_size;k++)
 198               fast[k]=x[k]*x[k];
 199           }
 200
 201
 202           if(sset->linkp==0 || firstlink==1){
 203
 204             compute_iir_freefall_limited(fast, input_size, &sss->iirS[i][j],
 205                                          smooth,smoothlimit);
 206
 207             memcpy(slow,fast,sizeof(slow));
 208             compute_iir_freefallonly1(slow, input_size, &sss->iirR[i][j],
 209                                       release);
 210
 211             //_analysis("fast3",i,fast,input_size,1,offset);
 212
 213             if(multiplier==sss->prevratio[i]){
 214
 215               for(k=0;k<input_size;k++)
 216                 fast[k]=fromdB_a((todB_a(slow[k])-todB_a(fast[k]))*.5*multiplier);
 217
 218             }else{
 219               float multiplier_add=(multiplier-sss->prevratio[i])/input_size;
 220               multiplier=sss->prevratio[i];
 221
 222               for(k=0;k<input_size;k++){
 223                 fast[k]=fromdB_a((todB_a(slow[k])-todB_a(fast[k]))*.5*multiplier);
 224                 multiplier+=multiplier_add;
 225               }
 226
 227             }
 228
 229             //_analysis("adj3",i,fast,input_size,1,offset);
 230
 231             if(sset->linkp && firstlink==1){
 232
 233               for(l=0;l<input_ch;l++){
 234                 if(l!=j){
 235                   memcpy(&sss->iirS[i][l],&sss->iirS[i][j],sizeof(iir_state));
 236                   memcpy(&sss->iirR[i][l],&sss->iirR[i][j],sizeof(iir_state));
 237                 }
 238               }
 239             }
 240
 241             firstlink++;
 242           }
 243
 244
 245           {
 246             float *x=sss->ss.lap[i][j];
 247             for(k=0;k<input_size;k++)
 248               if(fast[k]<1.)
 249                 x[k]*=fast[k];
 250           }
 251         }
 252       }
 253     }
 254
 255     sss->prevratio[i]=multiplier;
 256
 257   }
 258
 259 }
 260
 261 static void deverb_work_channel(void *vs){
 262   deverb_work_helper(vs,&deverb_channel_set);
 263 }
 264
 265 time_linkage *deverb_read_channel(time_linkage *in){
 266   int visible[input_ch];
 267   int active [input_ch];
 268   subband_window *w[input_ch];
 269   int i;
 270
 271   for(i=0;i<input_ch;i++){
 272     visible[i]=0;
 273     active[i]=deverb_channel_set.active[i];
 274     w[i]=&sw;
 275   }
 276
 277   return subband_read(in, &channel_state.ss, w, visible, active,
 278                       deverb_work_channel, &channel_state);
 279 }