function M = chol_ts(A)

%This program does Cholesky Decomposition of SPD matrices 

%mbrealvector(A);

[n,junk] = size(A);
   
M = zeros(n,n);
   
k = 1 ;

while k <= n 
   
     
   A(k,k) = sqrt( A(k,k) - sum(A(k,1:k-1).^2) ) ;
   
   A(1:k-1,k) = zeros(k-1,1);
   
   i = k + 1 ;
   
   while i <= n
      
      A(i,k) = (A(i,k) - (A(i,1:k-1)*A(k,1:k-1)'))/A(k,k); 
      
      i = i + 1; 
      
   end
   
   k = k + 1;
   
end

M = A;