how to plot a 2D covariance error ellipse?

Question

I have a location of landmark in 2D. According to Extended Kalman Filter EKF- SLAM, if the robot re-observes the same landmark, the covariance ellipse will shrink. I collected the necessary information and I would like to know how the covariance ellipse is drawn. The location of a landmark is $<\!x:30,y:60\!>$. Now for the first time the robot detects the location the following information is gathered.

$$ \mu_{x} = 28.8093 \\ \mu_{y} = 60.6267 \\ Cov(x,y) = \begin{bmatrix} 1.68165 & -0.793713 \\ -0.793713 & 0.388516 \\ \end{bmatrix} $$

I stored all the values in txt file. What is the formula for drawing the covariance?

Some samples from the experiment.

\mu_{x}     \mu_{y}     \sigma_{xx}  \sigma_{xy}    \sigma_{yx}     \sigma_{yy}  
---------------------------------------------------------------------------------
28.8093     60.6267     1.68165      -0.793713      -0.793713       0.388516
29.0079     60.5671     1.56697      -0.740083      -0.740083       0.358862
29.0511     60.5439     1.54802      -0.732739      -0.732739       0.353890
29.0148     60.5132     1.54433      -0.732171      -0.732171       0.352841
28.9692     60.4775     1.54340      -0.732399      -0.732399       0.352388
28.948      60.4577     1.54311      -0.732623      -0.732623       0.352052
28.9527     60.4621     1.54300      -0.732781      -0.732781       0.351782
28.9514     60.4602     1.54290      -0.732913      -0.732913       0.351591
28.9506     60.4596     1.54283      -0.733016      -0.733016       0.351445
28.9474     60.4539     1.54279      -0.733090      -0.733090       0.351320

---

Edit:

I have found this Matlab Code for drawing what I'm looking for but I don't understand the rule of Choleski method in the code.

NP = 16;
alpha  = 2*pi/NP*(0:NP);
circle = [cos(alpha);sin(alpha)];
ns = 3; 
x = [28.8093 ;60.626];
P = [1.68165 -0.793713;-0.793713 0.388516];
C = chol(P)'; %Choleski method <-????????????
ellip = ns*C*circle;
X = x(1)+ellip(1,:);
Y = x(2)+ellip(2,:);

The result is in the below picture which is exactly what I'm looking for but what is the rule of Choleski method in the code?

Answer 1

Given a mean vector and a covariance matrix, you could do something like the following:

library(MASS)
reps = mvrnorm(10000,c(28.8093,60.6267),Sigma=matrix(c(1.68165,-0.793713,-0.7937130,.388516),2))
dd = kde2d(reps[,1],reps[,2],n=200)
image(dd)
z = dd$z/sum(dd$z)
l = sort(z)[sum(cumsum(sort(z))<.05)]#this is the maximum at 5%
contour(x = dd$x,y=dd$y, z=z,level=l,add=T,lty=2)

This is a bit of a hack, there are probably better nonparametric approaches (i.e.: that don't assume multivariate normality. Also, the above assumes that the covariance is the covariance of the sample means, not of the measurements. If it is the latter, then what I plotted was just the covariance of the measurements, not a confidence interval around the mean.