double pinterp(int n,double*x,double*y,double z){
// Lagrange interpolating polynomial
	double s=0.0;
	for(int i=0; i<n; i++){
		double p=1.0;
		for(int k=0; k<n; k++) if(k!=i) p*=(z-x[k])/(x[i]-x[k]);
		s+=y[i]*p;
		}
	return s;
}