ftp.nice.ch/pub/next/unix/audio/cmix.s.tar.gz#/cmix/lpc/alpol.c

This is alpol.c in view mode; [Download] [Up]

#include "lpc.h"
alpol(sig, errn, rms1, rms2, c)
double *sig, *errn, *rms1, *rms2, *c;
{
        double a[POLEMAX][POLEMAX], v[POLEMAX], b[POLEMAX];
	double x[FRAMAX], y[FRAMAX];
        double *vp=v, *bp=b, *xp=x, *yp=y;
        double sum, sumx, sumy;
        int k1, i, l, k, limit, j;

        for (xp=x; xp-x < FRAME ;++xp,++sig) 
                *xp = *sig;
        k1 = NPOLE + 1;
        for (i=0; i < NPOLE ;++i)  {
                sum = (double) 0.0;
                for (k=NPOLE; k < FRAME ;++k)
                        sum += x[k-(i+1)] * x[k];
                v[i] = -sum;
                if (i != NPOLE - 1)  {
                        limit = NPOLE - (i+1);
                        for (l=0; l < limit ;++l)  {
                            sum += x[NPOLE-(i+1)-(l+1)]* x[NPOLE-(l+1)] - x[FRAME-(i+1)-(l+1)]* x[FRAME-(l+1)];
                            a[(i+1)+l][l] = a[l][(i+1)+l] = sum;
                        }
                }
        }
        sum = (double) 0.0;
        for (k=NPOLE; k < FRAME ;++k)
                sum += pow(x[k], (double) 2.0);
        sumy = sumx = sum;
        for (l=0; l < NPOLE ;++l)  {
                sum += pow(x[NPOLE-(l+1)], (double) 2.0) - pow(x[FRAME-(l+1)], (double) 2.0);
                a[l][l] = sum;
        }
        gauss(a, v, b);
        filtn(x, y, b);
        for (i=0; i < NPOLE ;++i)
                sumy = sumy - b[i]*v[i];
        *rms1 = sqrt(sumx/(double) (FRAME - k1 + 1) );
        *rms2 = sqrt(sumy/(double) (FRAME - k1 + 1) );
        *errn = pow(((*rms2)/(*rms1)), (double) 2.0);
        for (bp=b; bp-b < NPOLE ;++bp,++c)
                *c = *bp;
        return(0);
}

These are the contents of the former NiCE NeXT User Group NeXTSTEP/OpenStep software archive, currently hosted by Netfuture.ch.