1
2 """
3 covellipse package
4 Author: Joris De Ridder
5
6 Sometime people plot an "errorbox" in a 2D diagram while they should actually
7 plot an "error-ellipse". This package allows to draw such an error-ellipse
8 given the 2D covariance matrix of your measurement. The latter also allows
9 to take into account the correlation between the two quantitites.
10
11 Example:
12
13 Import the required packages:
14
15 >>> from numpy import *
16 >>> from numpy.random import multivariate_normal
17 >>> from matplotlib import pylab as p
18 >>> from covellipse import sigmaEllipse
19
20 For the sake of this example, we make some fake data:
21
22 >>> covMatrix = array([[2.0, -0.75*sqrt(4)*sqrt(2)],[-0.75*sqrt(4)*sqrt(2), 4.0]])
23 >>> m = multivariate_normal([0,0],covMatrix, 3000)
24
25 Compute the 2-sigma ellipse, centered around (0,0):
26
27 >>> x,y=sigmaEllipse(0,0,covMatrix,2)
28
29 Make a plot of the data as well as the ellipse:
30
31 >>> p.scatter(m[:,0], m[:,1], s=5)
32 >>> p.plot(x,y,linewidth=2)
33 """
34
35
36 import numpy as np
37 from math import atan,sqrt,pi
38
39
41
42 """
43 Returns x and y values of a complete ellipse
44
45 @param xc: x-coordinate of the center of the ellipse
46 @type xc: float
47 @param yc: y-coordinate of the center of the ellipse
48 @type yc: float
49 @param a: half the long axis of the ellipse
50 @type a: float
51 @param b: half the small axis of the ellipse
52 @type b: float
53 @param phi: angle between long axis and the x-axis (radians)
54 @type phi: float
55 @return: x,y: coordinate arrays of ellipse
56 @rtype: tuple of 2 arrays (dtype: float)
57 """
58
59 t = np.linspace(0,2*pi,100)
60 x = xc + a * np.cos(t) * np.cos(phi) - b * np.sin(t) * np.sin(phi)
61 y = yc + a * np.cos(t) * np.sin(phi) + b * np.sin(t) * np.cos(phi)
62
63 return x,y
64
65
66
67
69
70 """
71 Return x and y value of the |nSigma|-sigma ellipse
72 Given a covariance matrix, and the center of the ellipse (xc,yc)
73
74 @param xc: x-coordinate of the center of the ellipse
75 @type xc: float
76 @param yc: y-coordinate of the center of the ellipse
77 @type yc: float
78 @param covMatrix: 2D covariance matrix
79 @type covMatrix: 2x2 float nd-array
80 @param nSigma: the nSigma-contour of the ellipse will be returned
81 @type nSigma: float
82 @return: (x,y): the x and y values of the ellipse
83 @rtype: tuple of two float ndarrays
84 """
85
86 a = covMatrix[0,0]
87 b = covMatrix [0,1]
88 c = b
89 d = covMatrix[1,1]
90
91 eigenvalue1 = 0.5*(a+d+np.sqrt((a-d)**2+4*b*c))
92 eigenvalue2 = 0.5*(a+d-np.sqrt((a-d)**2+4*b*c))
93
94 semimajor = sqrt(eigenvalue1) * nSigma
95 semiminor = sqrt(eigenvalue2) * nSigma
96
97 if b == 0.0:
98 eigenvector1 = np.array([1, 0])
99 eigenvector2 = np.array([0, 1])
100 phi = 0.0
101 else:
102 eigenvector1 = np.array([b, eigenvalue1-a])
103 eigenvector2 = np.array([b, eigenvalue2-a])
104 phi = atan((eigenvalue1-a)/b)
105
106 x,y = ellipse(xc, yc, semimajor, semiminor, phi)
107
108 return x,y
109