function image_3d = interp_3d(welec,xelec,yelec,zelec,xs,ys,zs,p,q)
%	image = tstfit(welec,xelec,yelec,zelec,xs,ys,zs,p,q)
%	
n = size(xelec,2);
ms = size(xs,1);
ns = size(xs,2);
image_3d = zeros(ms,ns);
qm = zeros(ms,ns);
qm = q(1) +q(2)*xs +q(3)*ys+q(4)*xs.^2+q(5)*xs.*ys +q(6)*ys.^2 +q(7)*zs +q(8)*zs.*xs+q(9)*zs.*ys+q(10)*zs.^2;
sum = zeros(ms,ns);
s = zeros(ms,ns);
t = zeros(ms,ns);
r=zeros(ms,ns);
str = zeros(ms,ns);
for i=1:n
	s = xs - xelec(i); %*ones(ms,ns);
	t = ys - yelec(i); %*ones(ms,ns);
	r = zs - zelec(i); %*ones(ms,ns);
	str = s.^2+t.^2+r.^2;
	sum = sum + p(i)*((str+welec).^2).*log(str+welec);
end;
image_3d = sum+qm;

