float warppol(sig,past,d,c)
float sig,*past,d,*c;
{
	float temp1,temp2;
	int n;

	temp1 = past[NPOLEM1];
	past[NPOLEM1] = cq * outold - d * past[NPOLEM1];
	for(n=NPOLE-2; n>=0; n--) {
		temp2 = past[n];
		past[n] = d * (past[n+1] - past[n]) + temp1;
		temp1 = temp2;
		}
	for(n=0;n<NPOLE;n++)  sig += c[n] * past[n];
	outold = sig;
	return(sig);
}
float warpset(d,c)     
float d,*c;
{
	int m;
	float cl;

	for(m=1; m<NPOLE; m++)   c[m] += d * c[m-1];
	cl = 1./(1.-d * c[NPOLEM1]);
	cq = cl * (1. - d * d);
	return(cl);
}
warpinit()
{
	outold = 0;
}
bwarppol(sig,past,d,c,out,nvals)
float *sig,*past,d,*c,*out;
{
	float temp1,temp2;
	int i,n;
	for(i=0; i<nvals; i++) {
		temp1 = past[NPOLEM1];
		past[NPOLEM1] = cq * outold - d * past[NPOLEM1];
		for(n=NPOLE-1; n>=0; n--) {
			temp2 = past[n];
			past[n] = d * (past[n+1] - past[n]) + temp1;
			temp1 = temp2;
			}
		for(n=0;n<NPOLE;n++)   *sig += c[n] * past[n];
		*out++ = outold = *sig++;
	}
}