#include "/musr/H/ugens.h"
#include "/musr/H/sfheader.h"
#include <stdio.h>
#define RESIDAMP     0
#define RMSAMP       1
#define THRESH       2
#define PITCH	     3
			/* values in locs 0-3 of frame*/
#define NPOLE        32 
#define NPOLEM1      31
#define FRAMSIZE     (NPOLE+4)
#define FPREC        22
#define RECSIZE      (FPREC*FRAMSIZE)
#define BPREC        (RECSIZE*FLOAT)
#define BPFRAME      (FRAMSIZE*FLOAT)

int anal;
float cq,outold;
extern SFHEADER sfdesc[NFILES];
int anal;
char dataset_name[80];

lpcplay(p,n_args)
float *p;
{
	float amp,si,hn,phs,*f,srd2,magic,d,warpset();
	float c[FRAMSIZE],past[NPOLE*2],frames,frame1,frameno,ampmlt,errno;
	float ballpole(),alpvals[2048],buzvals[4096],pchval[2048];
	float thresh,randamp,randoff;
	float *cpoint;
	float x,transposition,newpch;
	int jcount,seed,i,nsamps,counter;
	float cps,tblvals[2],weight();
	int input, output, j,k;


	if(!n_args) {
		/*
		printf(" p[0]=starting time, p[1]=duration, p[2]=pitch,p[3]=frame1, p4=frame2 p[5]=amp, p[6]=thresh, p[7]=d,p[8]=inputskip\n");
		*/
		return;
	}
	for(i=0; i<NPOLE*2; i++) past[i] = 0;
	input = 0; output=1;

	srd2 = SR/2.;
	magic = 512./SR;
	jcount = phs = counter = 0; 
	randamp = .1;
	seed = .1;
	cpoint = c + 4;

/* p0=start,p1=dur,p2=8ve.pch,p3=frame1,p4=frame2,p5=amp,p6=thresh,p7=d,p8=inputskip*/
	nsamps = setnote(p[0],p[1],output);
	       	setnote(p[8],p[1],input);
	thresh = p[6];
	frames = p[4] - p[3] + 1.;
	frame1 = p[3];
	for(i=p[3]; i<= p[4]; i++) {
		d = i;
		getfr(d,c);
		pchval[i - (int)p[3]] = c[PITCH];
		if(c[PITCH] == 0) pchval[i - (int)p[3]]=55;
	}
	if(ABS(p[2]) < 1.)  transposition = octpch(p[2]);
	else	transposition = 
		octpch(p[2]) - octcps(weight(p[2],p[3],p[4],thresh));
printf("%f\n",transposition);
	tableset(p[1],(int)(p[4]-p[3]+1),tblvals);
	amp = p[5];
	d = p[7];
	f = (float *)floc(1);

	if(anal <=0 ) {
		printf("No open dataset!");
		closesf();
	}

	warpinit();
	for(i=nsamps; i>0; i -=counter) {
			frameno = ((float)(nsamps - i)/nsamps)*frames + frame1;
			if(getfr(frameno,c) == -1) break;
			errno = (c[2] > thresh) ? 0 : 1;  
			ampmlt = (errno) ? amp * c[0] : amp * c[0] * randamp;
			/*
			cps = tablei(nsamps-i,pchval,tblvals);
			newpch = cpsoct(octcps(cps)+transposition);
			if(newpch == 0) newpch=55;
			*/
			/* CRIPPLED THIS CODE SO NEVER USES PCH */
			newpch=112;
			si = newpch * magic;
			hn = (int)(srd2/newpch);
			counter = (float)(SR/newpch);
			counter = (counter > i) ? i : counter;
			if(counter <= 0) break;

			bgetin(buzvals,input,counter*sfchans(&sfdesc[input]));
			if(sfchans(&sfdesc[input]) == 2)
				for(j=0,k=0; j<counter; j++,k += 2) 
					buzvals[j] = buzvals[k]+buzvals[k+1];
			if(d) {
                                ampmlt *=  warpset(d,cpoint);
                                bwarppol(buzvals,past,d,cpoint,alpvals,counter);
			}
			else
				ballpole(buzvals,&jcount,NPOLE,
		               			past,cpoint,alpvals,counter);

			bwipeout(alpvals,output,counter);   
		}
	endnote(output);
}
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-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];
		*out++ = outold = *sig++;
	}
}
getfr(frameno,c)
float frameno,*c;
{
	int frame,i,j;
	static float array[RECSIZE];
	float fraction;
	static int oldframe = 0;
	static int endframe = 0;
	frame = (int)frameno;
	fraction = frameno - (float)frame;
	if(!((frame >= oldframe) && (frame < endframe))) {
		if(lseek(anal,((long)frame*(long)BPFRAME),0) == -1) {
			fprintf(stderr,"bad lseek on analysis file \n");
			return(-1);
		} 
		if(read(anal,(char *)array,BPREC) <= 0) {
			fprintf(stderr,"reached eof on analysis file \n");
			return(-1);
		}
		oldframe = frame;
		endframe = oldframe + FPREC - 1;
		}

	for(i=(frame-oldframe)*FRAMSIZE,j=0; j<FRAMSIZE; i++,j++)  
		*(c+j) = *(array+i) + fraction * (*(array+i+FRAMSIZE) 
						- *(array+i));
	return(0);
}
float table(nsample,array,tab)

long nsample;
float *array,*tab;

{
	register loc = ((float)(nsample)/(*tab)) * *(tab+1);
	if(loc >= *(tab+1)) loc = *(tab+1);
	return(*(array + loc));
}
float weight(newpch,frame1,frame2,thresh)
float newpch,frame1,frame2;
{
	float c[FRAMSIZE];
	int i;
	float xweight,sum;
	xweight = sum = 0;
	for(i=(int)frame1; i<(int)frame2; i++) {
		getfr((float)i,c);
		if(c[THRESH] > thresh) continue;
		xweight += c[RMSAMP];
		sum += (c[PITCH] * c[RMSAMP]);
	}
	return(sum/xweight);
}
dataset(p,n_args)
/* p1=dataset name, p2=npoles */
float *p;
{
	char *name;
	int i;

	i=(int)p[0];
	name=(char *)i;
	if(strcmp(name,dataset_name)== 0) {
		printf("\n%s is already open.\n",name);
		return;
		}
	strcpy(dataset_name,name);
	if((anal = open(name,0)) <= 0) {
		printf("Can't open %s\n",name);
		closesf();
		}
	printf("\nOpened dataset %s.\n",name);
}