Commit 59270848 authored by Monty Montgomery's avatar Monty Montgomery
Browse files

Must... get... changes... under... source... control...

Subversion is not scary, and there are no monsters under my bed.

Monty



git-svn-id: https://svn.xiph.org/trunk/postfish@6506 0101bb08-14d6-0310-b084-bc0e0c8e3800
parent a6e4c7e4
......@@ -2,7 +2,8 @@
# Fuck the horse it rode in on
# and Fuck its little dog Libtool too
CC=gcc
ADD_DEF= -DUGLY_IEEE754_FLOAT32_HACK=1
CC=gcc
LD=gcc
INSTALL=install
PREFIX=/usr/local
......@@ -10,34 +11,26 @@ BINDIR=$PREFIX/bin
ETCDIR=/etc
MANDIR=$PREFIX/man
# is this a platform that uses IEEE 754/854 32 bit floats? The
# following is good for a speedup on most of these systems, otherwise
# comment it out. Using this define on a system where a 'float' is
# *not* an IEEE 32 bit float will destroy, destroy, destroy the audio.
IEEE=-DNASTY_IEEE_FLOAT32_HACK_IS_FASTER_THAN_LOG=1
SRC = main.c mainpanel.c multibar.c readout.c input.c output.c clippanel.c \
declip.c reconstruct.c multicompand.c windowbutton.c subpanel.c \
feedback.c freq.c eq.c eqpanel.c compandpanel.c subband.c lpc.c \
bessel.c
bessel.c suppresspanel.c suppress.c singlecomp.c singlepanel.c \
limit.c limitpanel.c
OBJ = main.o mainpanel.o multibar.o readout.o input.o output.o clippanel.o \
declip.o reconstruct.o multicompand.o windowbutton.o subpanel.o \
feedback.o freq.o eq.o eqpanel.o compandpanel.o subband.o lpc.o \
bessel.o
bessel.o suppresspanel.o suppress.o singlecomp.o singlepanel.o \
limit.o limitpanel.o
GCF = `pkg-config --cflags gtk+-2.0` -DG_DISABLE_DEPRECATED -DGDK_DISABLE_DEPRECATED -DGTK_DISABLE_DEPRECATED -DGDK_PIXBUF_DISABLE_DEPRECATED
all:
$(MAKE) target CFLAGS="-W -O3 -ffast-math $(GCF) $(IEEE)"
$(MAKE) target CFLAGS="-O3 -ffast-math -fomit-frame-pointer $(GCF) $(ADD_DEF)"
debug:
$(MAKE) target CFLAGS="-g -W -D__NO_MATH_INLINES $(GCF) $(IEEE)"
$(MAKE) target CFLAGS="-g -Wall -W -Wno-unused-parameter -D__NO_MATH_INLINES $(GCF) $(ADD_DEF)"
profile:
$(MAKE) target CFLAGS="-W -pg -g -O3 -ffast-math $(GCF) $(IEEE)" LIBS="-lgprof-helper"
$(MAKE) target CFLAGS="-pg -g -O3 -ffast-math $(GCF) $(ADD_DEF)" LIBS="-lgprof-helper"
clean:
rm -f $(OBJ) *.d gmon.out
......
......@@ -23,20 +23,52 @@
/* code derived directly from mkfilter by the late A.J. Fisher,
University of York <fisher@minster.york.ac.uk> September 1992; this
is only the minimum code needed to build an arbitrary 2nd order
Bessel filter */
is only the minimum code needed to build an arbitrary Bessel filter */
#include "postfish.h"
#define TWOPI (2.0 * M_PIl)
#define EPS 1e-10
#define MAXORDER 2
#define MAXPZ 4
#include "bessel.h"
typedef struct {
double re, im;
} complex;
static complex bessel_poles[] = {
{ -1.00000000000e+00, 0.00000000000e+00},
{ -1.10160133059e+00, 6.36009824757e-01},
{ -1.32267579991e+00, 0.00000000000e+00},
{ -1.04740916101e+00, 9.99264436281e-01},
{ -1.37006783055e+00, 4.10249717494e-01},
{ -9.95208764350e-01, 1.25710573945e+00},
{ -1.50231627145e+00, 0.00000000000e+00},
{ -1.38087732586e+00, 7.17909587627e-01},
{ -9.57676548563e-01, 1.47112432073e+00},
{ -1.57149040362e+00, 3.20896374221e-01},
{ -1.38185809760e+00, 9.71471890712e-01},
{ -9.30656522947e-01, 1.66186326894e+00},
{ -1.68436817927e+00, 0.00000000000e+00},
{ -1.61203876622e+00, 5.89244506931e-01},
{ -1.37890321680e+00, 1.19156677780e+00},
{ -9.09867780623e-01, 1.83645135304e+00},
{ -1.75740840040e+00, 2.72867575103e-01},
{ -1.63693941813e+00, 8.22795625139e-01},
{ -1.37384121764e+00, 1.38835657588e+00},
{ -8.92869718847e-01, 1.99832584364e+00},
{ -1.85660050123e+00, 0.00000000000e+00},
{ -1.80717053496e+00, 5.12383730575e-01},
{ -1.65239648458e+00, 1.03138956698e+00},
{ -1.36758830979e+00, 1.56773371224e+00},
{ -8.78399276161e-01, 2.14980052431e+00},
{ -1.92761969145e+00, 2.41623471082e-01},
{ -1.84219624443e+00, 7.27257597722e-01},
{ -1.66181024140e+00, 1.22110021857e+00},
{ -1.36069227838e+00, 1.73350574267e+00},
{ -8.65756901707e-01, 2.29260483098e+00},
};
#define TWOPI (2.0 * M_PIl)
#define EPS 1e-10
#define MAXPZ 8
typedef struct {
complex poles[MAXPZ], zeros[MAXPZ];
int numpoles, numzeros;
......@@ -66,6 +98,11 @@ static complex csub(complex z1, complex z2){
return z1;
}
static complex cconj(complex z){
z.im = -z.im;
return z;
}
static complex eval(complex coeffs[], int npz, complex z){
complex sum = (complex){0.0,0.0};
int i;
......@@ -105,8 +142,8 @@ static void expand(complex pz[], int npz, complex coeffs[]){
}
}
double mkbessel_2(double raw_alpha,double *ycoeff0,double *ycoeff1){
int i;
double mkbessel(double raw_alpha,int order,double *ycoeff){
int i,p= (order*order)/4;
pzrep splane, zplane;
complex topcoeffs[MAXPZ+1], botcoeffs[MAXPZ+1];
double warped_alpha;
......@@ -115,8 +152,11 @@ double mkbessel_2(double raw_alpha,double *ycoeff0,double *ycoeff1){
memset(&splane,0,sizeof(splane));
memset(&zplane,0,sizeof(zplane));
splane.poles[splane.numpoles++] = (complex){ -1.10160133059e+00, 6.36009824757e-01};
splane.poles[splane.numpoles++] = (complex){ -1.10160133059e+00, -6.36009824757e-01};
if (order & 1) splane.poles[splane.numpoles++] = bessel_poles[p++];
for (i = 0; i < order/2; i++){
splane.poles[splane.numpoles++] = bessel_poles[p];
splane.poles[splane.numpoles++] = cconj(bessel_poles[p++]);
}
warped_alpha = tan(M_PIl * raw_alpha) / M_PIl;
for (i = 0; i < splane.numpoles; i++){
......@@ -137,11 +177,296 @@ double mkbessel_2(double raw_alpha,double *ycoeff0,double *ycoeff1){
expand(zplane.poles, zplane.numpoles, botcoeffs);
dc_gain = evaluate(topcoeffs, zplane.numzeros, botcoeffs, zplane.numpoles, (complex){1.0,0.0});
*ycoeff0 = -(botcoeffs[0].re / botcoeffs[zplane.numpoles].re);
*ycoeff1 = -(botcoeffs[1].re / botcoeffs[zplane.numpoles].re);
for(i=0;i<order;i++)
ycoeff[order-i-1] = -(botcoeffs[i].re / botcoeffs[zplane.numpoles].re);
return hypot(dc_gain.re,dc_gain.im);
}
/* assymetrical attack/decay filter computation */
/* this one is designed for fast attack, slow decay */
void compute_iir_fast_attack2(float *x, int n, iir_state *is,
iir_filter *attack, iir_filter *decay){
double a_c0=attack->c[0],d_c0=decay->c[0];
double a_c1=attack->c[1],d_c1=decay->c[1];
double a_g=attack->g, d_g=decay->g;
double x0=is->x[0],x1=is->x[1];
double y0=is->y[0],y1=is->y[1];
int state=is->state;
int i=0;
if(x[0]>y0)state=0;
while(i<n){
if(state==0){
/* attack case */
while(i<n){
double ya= (x[i]+x0*2.+x1)/a_g + y0*a_c0+y1*a_c1;
if(ya<y0){
state=1;
break;
}
x1=x0;x0=x[i];
y1=y0;x[i]=y0=ya;
i++;
}
}
if(state==1){
/* decay case */
if(y1<y0){
/* decay fixup needed because we're in discontinuous time */
y1=y0;
}
while(1){
double yd = (x[i]+x0*2.+x1)/d_g + y0*d_c0+y1*d_c1;
x1=x0;x0=x[i];
y1=y0;x[i]=y0=yd;
i++;
if(i>=n)break;
if(x[i]>y0){
state=0;
break;
}
}
}
}
is->x[0]=x0;is->x[1]=x1;
is->y[0]=y0;is->y[1]=y1;
is->state=state;
}
/* allow decay to proceed in freefall */
void compute_iir_freefall1(float *x, int n, iir_state *is,
iir_filter *decay){
double d_c0=decay->c[0];
double x0=is->x[0];
double y0=is->y[0];
int i=0;
while(i<n){
double yd;
yd = y0*d_c0;
if(x[i]>yd)yd=x[i];
x0=x[i];
x[i]=y0=yd;
i++;
}
is->x[0]=x0;
is->y[0]=y0;
}
void compute_iir_freefall2(float *x, int n, iir_state *is,
iir_filter *decay){
double d_c0=decay->c[0];
double d_c1=decay->c[1];
double x0=is->x[0];
double x1=is->x[1];
double y0=is->y[0];
double y1=is->y[1];
int i=0;
while(i<n){
double yd;
if(y1<y0)y1=y0; // slope fixup
yd = y0*d_c0+y1*d_c1;
if(x[i]>yd)yd=x[i];
x1=x0;x0=x[i];
y1=y0;x[i]=y0=yd;
i++;
}
is->x[0]=x0;
is->x[1]=x1;
is->y[0]=y0;
is->y[1]=y1;
}
void compute_iir_freefall3(float *x, int n, iir_state *is,
iir_filter *decay){
double d_c0=decay->c[0];
double d_c1=decay->c[1];
double d_c2=decay->c[2];
double x0=is->x[0],y0=is->y[0];
double x1=is->x[1],y1=is->y[1];
double x2=is->x[2],y2=is->y[2];
int i=0;
while(i<n){
double yd;
if(y1<y0)y1=y0; // slope fixup
if(y2<y1)y2=y1; // slope fixup
yd = y0*d_c0+y1*d_c1+y2*d_c2;
if(x[i]>yd)yd=x[i];
x2=x1;x1=x0;x0=x[i];
y2=y1;y1=y0;x[i]=y0=yd;
i++;
}
is->x[0]=x0;is->y[0]=y0;
is->x[1]=x1;is->y[1]=y1;
is->x[2]=x2;is->y[2]=y2;
}
void compute_iir_freefall4(float *x, int n, iir_state *is,
iir_filter *decay){
double d_c0=decay->c[0];
double d_c1=decay->c[1];
double d_c2=decay->c[2];
double d_c3=decay->c[3];
double x0=is->x[0],y0=is->y[0];
double x1=is->x[1],y1=is->y[1];
double x2=is->x[2],y2=is->y[2];
double x3=is->x[3],y3=is->y[3];
int i=0;
while(i<n){
double yd;
if(y1<y0)y1=y0; // slope fixup
if(y2<y1)y2=y1; // slope fixup
if(y3<y2)y3=y2; // slope fixup
yd = y0*d_c0+y1*d_c1+y2*d_c2+y3*d_c3;
if(x[i]>yd)yd=x[i];
x3=x2;x2=x1;x1=x0;x0=x[i];
y3=y2;y2=y1;y1=y0;x[i]=y0=yd;
i++;
}
is->x[0]=x0;is->y[0]=y0;
is->x[1]=x1;is->y[1]=y1;
is->x[2]=x2;is->y[2]=y2;
is->x[3]=x3;is->y[3]=y3;
}
/* symmetric filter computation */
void compute_iir_symmetric2(float *x, int n, iir_state *is,
iir_filter *filter){
double c0=filter->c[0];
double c1=filter->c[1];
double g=filter->g;
double x0=is->x[0];
double x1=is->x[1];
double y0=is->y[0];
double y1=is->y[1];
int i=0;
while(i<n){
double yd= (x[i]+x0*2.+x1)/g + y0*c0+y1*c1;
x1=x0;x0=x[i];
y1=y0;x[i]=y0=yd;
i++;
}
is->x[0]=x0;
is->x[1]=x1;
is->y[0]=y0;
is->y[1]=y1;
}
void compute_iir_symmetric3(float *x, int n, iir_state *is,
iir_filter *filter){
double c0=filter->c[0];
double c1=filter->c[1];
double c2=filter->c[2];
double g=filter->g;
double x0=is->x[0],y0=is->y[0];
double x1=is->x[1],y1=is->y[1];
double x2=is->x[2],y2=is->y[2];
int i=0;
while(i<n){
double yd= (x[i]+(x0+x1)*3.+x2)/g + y0*c0+y1*c1+y2*c2;
x2=x1;x1=x0;x0=x[i];
y2=y1;y1=y0;x[i]=y0=yd;
i++;
}
is->x[0]=x0;is->y[0]=y0;
is->x[1]=x1;is->y[1]=y1;
is->x[2]=x2;is->y[2]=y2;
}
void compute_iir_symmetric4(float *x, int n, iir_state *is,
iir_filter *filter){
double c0=filter->c[0];
double c1=filter->c[1];
double c2=filter->c[2];
double c3=filter->c[3];
double g=filter->g;
double x0=is->x[0],y0=is->y[0];
double x1=is->x[1],y1=is->y[1];
double x2=is->x[2],y2=is->y[2];
double x3=is->x[3],y3=is->y[3];
int i=0;
while(i<n){
double yd= (x[i]+(x0+x2)*4.+x1*6.+x3)/g +
y0*c0+y1*c1+y2*c2+y3*c3;
x3=x2;x2=x1;x1=x0;x0=x[i];
y3=y2;y2=y1;y1=y0;x[i]=y0=yd;
i++;
}
is->x[0]=x0;is->y[0]=y0;
is->x[1]=x1;is->y[1]=y1;
is->x[2]=x2;is->y[2]=y2;
is->x[3]=x3;is->y[3]=y3;
}
/* filter decision wrapper */
void compute_iir2(float *x, int n, iir_state *is,
iir_filter *attack, iir_filter *decay){
if (attack->alpha > decay->alpha){
/* fast attack, slow decay */
compute_iir_fast_attack2(x, n, is, attack, decay);
}else{
compute_iir_symmetric2(x, n, is, attack);
}
}
......@@ -21,4 +21,71 @@
*
*/
extern double mkbessel_2(double raw_alpha,double *ycoeff0,double *ycoeff1);
#include "postfish.h"
extern int input_rate;
#define MAXORDER 4
typedef struct {
double c[MAXORDER];
double g;
int order;
float alpha;
float Hz;
float ms;
} iir_filter;
static inline long impulse_ahead2(float alpha){
return rint(.13f/alpha);
}
static inline long impulse_ahead3(float alpha){
return rint(.22f/alpha);
}
static inline long impulse_ahead4(float alpha){
return rint(.32f/alpha);
}
static inline long step_ahead(float alpha){
return rint(.6f/alpha);
}
static inline float step_freq(long ahead){
return input_rate*.6f/ahead;
}
static inline float impulse_freq2(long ahead){
return input_rate*.13f/ahead;
}
static inline float impulse_freq3(long ahead){
return input_rate*.22f/ahead;
}
static inline float impulse_freq4(long ahead){
return input_rate*.32f/ahead;
}
typedef struct {
double x[MAXORDER];
double y[MAXORDER];
int state;
} iir_state;
extern double mkbessel(double raw_alpha,int order,double *ycoeff);
extern void compute_iir_fast_attack2(float *x, int n, iir_state *is,
iir_filter *attack, iir_filter *decay);
extern void compute_iir_symmetric2(float *x, int n, iir_state *is,
iir_filter *filter);
extern void compute_iir_symmetric3(float *x, int n, iir_state *is,
iir_filter *filter);
extern void compute_iir_symmetric4(float *x, int n, iir_state *is,
iir_filter *filter);
extern void compute_iir2(float *x, int n, iir_state *is,
iir_filter *attack, iir_filter *decay);
extern void compute_iir_freefall1(float *x, int n, iir_state *is,
iir_filter *decay);
extern void compute_iir_freefall2(float *x, int n, iir_state *is,
iir_filter *decay);
extern void compute_iir_freefall3(float *x, int n, iir_state *is,
iir_filter *decay);
extern void compute_iir_freefall4(float *x, int n, iir_state *is,
iir_filter *decay);
......@@ -59,12 +59,13 @@ static int bank_active=2;
static cbar bars[multicomp_freqs_max+1];
static int inactive_updatep=1;
static void under_compand_change(GtkWidget *w,gpointer in){
static void compand_change(GtkWidget *w,Readout *r,sig_atomic_t *var){
char buffer[80];
Readout *r=(Readout *)in;
float val=1./multibar_get_value(MULTIBAR(w),0);
if(val>=10){
if(val==1.){
sprintf(buffer," off");
}else if(val>=10){
sprintf(buffer,"%4.1f:1",val);
}else if(val>=1){
sprintf(buffer,"%4.2f:1",val);
......@@ -76,197 +77,73 @@ static void under_compand_change(GtkWidget *w,gpointer in){
readout_set(r,buffer);
c.under_ratio=rint(val*1000.);
*var=rint(val*1000.);
}
static void under_compand_change(GtkWidget *w,gpointer in){
compand_change(w,(Readout *)in,&c.under_ratio);
}
static void over_compand_change(GtkWidget *w,gpointer in){
char buffer[80];
Readout *r=(Readout *)in;
float val=1./multibar_get_value(MULTIBAR(w),0);
if(val>=10){
sprintf(buffer,"%4.1f:1",val);
}else if(val>=1){
sprintf(buffer,"%4.2f:1",val);
}else if(val>.10001){
sprintf(buffer,"1:%4.2f",1./val);
}else{
sprintf(buffer,"1:%4.1f",1./val);
}
compand_change(w,(Readout *)in,&c.over_ratio);
}
readout_set(r,buffer);
c.over_ratio=rint(val*1000.);
static void base_compand_change(GtkWidget *w,gpointer in){
compand_change(w,(Readout *)in,&c.base_ratio);
}
static void suppress_compand_change(GtkWidget *w,gpointer in){
static void timing_display(GtkWidget *w,Readout *r,float v){
char buffer[80];
Readout *r=(Readout *)in;
float val=1./multibar_get_value(MULTIBAR(w),0);
if(val>.10001){
sprintf(buffer,"1:%4.2f",1./val);
if(v<10){
sprintf(buffer,"%4.2fms",v);
}else if(v<100){
sprintf(buffer,"%4.1fms",v);
}else if (v<1000){
sprintf(buffer,"%4.0fms",v);
}else if (v<10000){
sprintf(buffer,"%4.2fs",v/1000.);
}else{
sprintf(buffer,"1:%4.1f",1./val);
sprintf(buffer,"%4.1fs",v/1000.);
}
readout_set(r,buffer);
c.suppress_ratio=rint(val*1000.);
}
static void under_limit_change(GtkWidget *w,gpointer in){
char buffer[80];
Readout *r=(Readout *)in;
float val=140.-multibar_get_value(MULTIBAR(w),0);
sprintf(buffer,"%3ddB",(int)rint(val));
readout_set(r,buffer);
c.under_limit=rint(val*10.);