import java.io.*;
import javax.swing.*;

class SCO15
{
public static double[] savgol(double c[], int nl, int nr, int ld, int m)
{
	int np=c.length;
	int imj,ipj,j,k,kk,mm;
	int indx[]=new int[m+1];
	double a[][]=new double[m+1][m+1];
	double b[]=new double[m+1];
	double d,fac,sum;
	if (np < nl+nr+1 || nl < 0 || nr < 0 || ld > m || nl+nr < m)
	System.out.println("bad args in savgol");
	for (ipj=0;ipj<=(m << 1);ipj++) {
		sum=((ipj>0) ? 0.0 : 1.0);
		for (k=1;k<=nr;k++) sum += Math.pow((double)k,(double)ipj);
		for (k=1;k<=nl;k++) sum += Math.pow((double)-k,(double)ipj);
		mm=Math.min(ipj,2*m-ipj);
		for (imj = -mm;imj<=mm;imj+=2) a[1+(ipj+imj)/2][1+(ipj-imj)/2]=sum;
	}
	d=0;
	ludcmp(a,d,indx);
	for (j=1;j<=m+1;j++) b[j]=0.0;
	b[ld+1]=1.0;
	lubksb(a,b,indx);
	for (kk=1;kk<=np;kk++) c[kk]=0.0;
	for (k = -nl;k<=nr;k++) {
		sum=b[1];
		fac=1.0;
		for (mm=1;mm<=m;mm++) sum += b[mm+1]*(fac *= k);
		kk=((np-k) % np)+1;
		c[kk]=sum;
	}
return c;
}

public static void lubksb(double [][]a, double b[],int indx[] )
{
	int n=a.length;
	int i,ii=0,ip,j;
	double sum;

	for (i=1;i<=n;i++) {
		ip=indx[i];
		sum=b[ip];
		b[ip]=b[i];
		if (ii>0)
			for (j=ii;j<=i-1;j++) sum -= a[i][j]*b[j];
		else if (sum>0) ii=i;
		b[i]=sum;
	}
	for (i=n;i>=1;i--) {
		sum=b[i];
		for (j=i+1;j<=n;j++) sum -= a[i][j]*b[j];
		b[i]=sum/a[i][i];
	}
}

public static void ludcmp(double a[][],double d, int indx[] )
{
	int n=a.length;
	int i,imax,j,k;
	imax=0;
	double big,dum,sum,temp;
	double vv[]=new double[n];
    double TINY=1e-10;
	d=1.0;
	for (i=1;i<=n;i++) {
		big=0.0;
		for (j=1;j<=n;j++)
			if ((temp=Math.abs(a[i][j])) > big) big=temp;
		if (big == 0.0) System.out.println("Singular matrix in routine ludcmp");
		vv[i]=1.0/big;
	}
	for (j=1;j<=n;j++) {
		for (i=1;i<j;i++) {
			sum=a[i][j];
			for (k=1;k<i;k++) sum -= a[i][k]*a[k][j];
			a[i][j]=sum;
		}
		big=0.0;
		for (i=j;i<=n;i++) {
			sum=a[i][j];
			for (k=1;k<j;k++)
				sum -= a[i][k]*a[k][j];
			a[i][j]=sum;
			if ( (dum=vv[i]*Math.abs(sum)) >= big) {
				big=dum;
				imax=i;
			}
		}
		if (j != imax) {
			for (k=1;k<=n;k++) {
				dum=a[imax][k];
				a[imax][k]=a[j][k];
				a[j][k]=dum;
			}
			d = -d;
			vv[imax]=vv[j];
		}
		indx[j]=imax;
		if (a[j][j] == 0.0) a[j][j]=TINY;
		if (j != n) {
			dum=1.0/(a[j][j]);
			for (i=j+1;i<=n;i++) a[i][j] *= dum;
		}
	}
}
public static void main(String arg[])
{
	int i,j,m,nl,np,nr;
	int NMAX=1000;
	int NTEST=6;
	double sum;
	int mtest[]={0,2,2,2,2,4,4};
	int nltest[]={0,2,3,4,5,4,5};
	int nrtest[]={0,2,1,0,5,4,5};
	String ans[]={"",
"                      -0.086  0.343  0.486  0.343 -0.086",
"               -0.143  0.171  0.343  0.371  0.257",
"         0.086 -0.143 -0.086  0.257  0.886",
" -0.084  0.021  0.103  0.161  0.196  0.207  0.196  0.161  0.103  0.021 -0.084",
"         0.035 -0.128  0.070  0.315  0.417  0.315  0.070 -0.128  0.035",
"  0.042 -0.105 -0.023  0.140  0.280  0.333  0.280  0.140 -0.023 -0.105  0.042"};

	double c[]=new double[NMAX];
	System.out.println("M nl nr\n");
	System.out.println("\t\t\tSample Savitzky-Golay Coefficients\n");
	for (i=1;i<=NTEST;i++) {
		m=mtest[i];
		nl=nltest[i];
		nr=nrtest[i];
		np=nl+nr+1;
		savgol(c,nl,nr,0,m);
		for (sum=0.0,j=1;j<=np;j++) sum += c[j];
		System.out.println("m="+m+"nl="+nl+"nr="+nr);
		for (j=nl;j<5;j++) System.out.println(" ");
		for (j=nl+1;j>=1;j--) System.out.println(c[j]);
		for (j=0;j<nr;j++) System.out.println(c[np-j]);
		System.out.println(" ");
		System.out.println("Sum = "+sum);
		System.out.println("Compare answer:"+ans[i]);
	}
}
}
