function pp=VVP_basis_dal(x,y,z,x0,y0,z0,t,p)

ae=6378.1e03;
r=sqrt(x.^2+y.^2+z.^2);
num=length(z);
zz=z0*ones(size(z));
xx=x0*ones(size(x));
yy=y0*ones(size(y));
thetae=t+atan((r.*cos(t))./(ae+r.*sin(t)));

iproj=cos(thetae).*sin(p);
jproj=cos(thetae).*cos(p);
kproj=sin(thetae);
xe=r.*iproj;
ye=r.*jproj;
ze=r.*kproj-((r.*cos(thetae)).^2./(2*ae));

pp(:,1)=iproj;
pp(:,2)=jproj;
pp(:,3)=kproj;
pp(:,4)=iproj.*(xe-xx);
pp(:,5)=iproj.*(ye-yy);
pp(:,6)=jproj.*(xe-xx);
pp(:,7)=jproj.*(ye-yy);
% pp(:,8)=kproj.*(xe-xx);
% pp(:,9)=kproj.*(ye-yy);

return
end