How to interpolate with more interpolation points?

Hi everyone, below this the script of interpolation. By using this script, I managed to interpolate if there is only one point / location where I want to estimate the data. However, I have a lot of points / locations to be estimated. I tried to modified the example script below, and the outcome for the first and second location are the same. Can anyone know what's the problem with the modified script?

clc;clear all; close all;

%% Collected Data
X = 1e3.*[40; -45; 45];
Y = 1e3.*[60;-30;-70];
Z = [24.6; 26.4; 25.2];

%% Location for data estimation
XI = 1e3*[10;11;20];
YI = 1e3*[15;16;19];

%% Covariance Data Matrix
for k=1:3
    for m=1:3
        dkm(k,m)=sqrt(abs(X(k)-X(m))^2+abs(Y(k)-Y(m))^2);
    end
end

Ckm=(1-0.01)*dkm

for k=1:3
    for j=1:length(XI)
    dki(k,length(XI))=sqrt(abs(X(k)-XI(j)).^2+abs(Y(k)-YI(j)).^2);
    end
end

Cki=(1-0.01)*dki

%%
w=inv(Ckm)*Cki

% BE CAREFUL: BAD ESTIMATION if means is not substracted
ZI=w'*Z

% GOOD ESTIMATION if the mean is removed
Zprim = Z-mean(Z);
ZIprim = w'*Zprim;
ZIgood = mean(Z)+ZIprim

Append2 = [XI YI ZI ZIgood];