function df=get_K(x,y,z,x_vis)

%precalculate all kernel matrices needed in algorithm

[x1, x2, y1, y2, ~, ~]=split(x,y,z,.5);

p=x(:,1);
t=x(:,2);
s=x(:,3);

vp=median_inter(p);
vt=median_inter(t);
vs=median_inter(s);

p1=x1(:,1);
t1=x1(:,2);
s1=x1(:,3);

p2=x2(:,1);
t2=x2(:,2);
s2=x2(:,3);

ptest=x_vis(:,1);
ttest=x_vis(:,2);
stest=x_vis(:,3);

df.y1=y1;
df.y2=y2;
df.y=y;

df.K11=get_K_matrix(p1,p1,vp).*get_K_matrix(t1,t1,vt).*get_K_matrix(s1,s1,vs);
df.K22=get_K_matrix(p2,p2,vp).*get_K_matrix(t2,t2,vt).*get_K_matrix(s2,s2,vs);
df.K12=get_K_matrix(p1,p2,vp).*get_K_matrix(t1,t2,vt).*get_K_matrix(s1,s2,vs);
df.Kaa=get_K_matrix(p,p,vp).*get_K_matrix(t,t,vt).*get_K_matrix(s,s,vs);
df.Kta=get_K_matrix(ptest,p,vp).*get_K_matrix(ttest,t,vt).*get_K_matrix(stest,s,vs);

%df.K11=get_K_matrix(p1,p1,vp)+get_K_matrix(t1,t1,vt)+get_K_matrix(s1,s1,vs);
%df.K22=get_K_matrix(p2,p2,vp)+get_K_matrix(t2,t2,vt)+get_K_matrix(s2,s2,vs);
%df.K12=get_K_matrix(p1,p2,vp)+get_K_matrix(t1,t2,vt)+get_K_matrix(s1,s2,vs);
%df.Kaa=get_K_matrix(p,p,vp)+get_K_matrix(t,t,vt)+get_K_matrix(s,s,vs);
%df.Kta=get_K_matrix(ptest,p,vp)+get_K_matrix(ttest,t,vt)+get_K_matrix(stest,s,vs);

end

