lib/psy.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-2000             *
   9  * by Monty <monty@xiph.org> and The XIPHOPHORUS Company            *
  10  * http://www.xiph.org/                                             *
  11  *                                                                  *
  12  ********************************************************************
  13
  14  function: psychoacoustics not including preecho
  15  last mod: $Id: psy.c,v 1.17 2000/03/10 13:21:18 xiphmont Exp $
  16
  17  ********************************************************************/
  18
  19 #include <stdlib.h>
  20 #include <math.h>
  21 #include <string.h>
  22 #include <stdio.h>
  23 #include "vorbis/codec.h"
  24
  25 #include "psy.h"
  26 #include "lpc.h"
  27 #include "smallft.h"
  28 #include "scales.h"
  29 #include "misc.h"
  30
  31 /* Set up decibel threshhold slopes on a Bark frequency scale */
  32
  33 static void set_curve(double *ref,double *c,int n, double crate){
  34   int i,j=0;
  35
  36   for(i=0;i<MAX_BARK-1;i++){
  37     int endpos=rint(fromBARK(i+1)*2*n/crate);
  38     double base=ref[i];
  39     double delta=(ref[i+1]-base)/(endpos-j);
  40     for(;j<endpos && j<n;j++){
  41       c[j]=base;
  42       base+=delta;
  43     }
  44   }
  45 }
  46
  47 void _vp_psy_init(vorbis_look_psy *p,vorbis_info_psy *vi,int n,long rate){
  48   long i;
  49   memset(p,0,sizeof(vorbis_look_psy));
  50   p->maskthresh=malloc(n*sizeof(double));
  51   p->barknum=malloc(n*sizeof(double));
  52   p->vi=vi;
  53   p->n=n;
  54
  55   /* set up the lookups for a given blocksize and sample rate */
  56   /* Vorbis max sample rate is limited by 26 Bark (54kHz) */
  57   set_curve(vi->maskthresh, p->maskthresh, n,rate);
  58
  59   for(i=0;i<n;i++)
  60     p->barknum[i]=toBARK(rate/2.*i/n);
  61
  62 #ifdef ANALYSIS
  63   {
  64     int j;
  65     FILE *out;
  66     char buffer[80];
  67
  68     sprintf(buffer,"mask_threshhold_%d.m",n);
  69     out=fopen(buffer,"w+");
  70     for(j=0;j<n;j++)
  71       fprintf(out,"%g\n",p->maskthresh[j]);
  72     fclose(out);
  73   }
  74 #endif
  75
  76 }
  77
  78 void _vp_psy_clear(vorbis_look_psy *p){
  79   if(p){
  80     if(p->maskthresh)free(p->maskthresh);
  81     if(p->barknum)free(p->barknum);
  82     memset(p,0,sizeof(vorbis_look_psy));
  83   }
  84 }
  85
  86 /* Masking curve: linear rolloff on a Bark/dB scale, attenuated by
  87    maskthresh */
  88
  89 void _vp_mask_floor(vorbis_look_psy *p,double *f, double *floor,int attp){
  90   int n=p->n;
  91   double hroll=p->vi->hrolldB;
  92   double lroll=p->vi->lrolldB;
  93   double curmask=todB(f[0])+(attp?p->maskthresh[0]:0);
  94   double curoc=0.;
  95   long i;
  96
  97   /* run mask forward then backward */
  98   for(i=0;i<n;i++){
  99     double newmask=todB(f[i])+(attp?p->maskthresh[i]:0);
 100     double newoc=p->barknum[i];
 101     double roll=curmask-(newoc-curoc)*hroll;
 102     double troll;
 103     if(newmask>roll){
 104       roll=curmask=newmask;
 105       curoc=newoc;
 106     }
 107     troll=fromdB(roll);
 108     if(floor[i]<troll)floor[i]=troll;
 109   }
 110
 111   curmask=todB(f[n-1])+(attp?p->maskthresh[n-1]:0);
 112   curoc=p->barknum[n-1];
 113   for(i=n-1;i>=0;i--){
 114     double newmask=todB(f[i])+(attp?p->maskthresh[i]:0);
 115     double newoc=p->barknum[i];
 116     double roll=curmask-(curoc-newoc)*lroll;
 117     double troll;
 118     if(newmask>roll){
 119       roll=curmask=newmask;
 120       curoc=newoc;
 121     }
 122     troll=fromdB(roll);
 123     if(floor[i]<troll)floor[i]=troll;
 124   }
 125 }
 126
 127 /* s must be padded at the end with m-1 zeroes */
 128 static void time_convolve(double *s,double *r,int n,int m){
 129   int i;
 130
 131   for(i=0;i<n;i++){
 132     int j;
 133     double acc=0;
 134
 135     for(j=0;j<m;j++)
 136       acc+=s[i+j]*r[m-j-1];
 137
 138     s[i]=acc;
 139   }
 140 }
 141
 142 void _vi_psy_free(vorbis_info_psy *i){
 143   if(i){
 144     memset(i,0,sizeof(vorbis_info_psy));
 145     free(i);
 146   }
 147 }