Rotating a gaussian function - matlab

Question

Good day,

I have code that stretches and rotates a gaussian 2d pdf as such:

mu = [0 0];

Sigma = [1 0; 0 1];
Scale = [3 0; 0 1];
Theta = 10;
Rotate = [cosd(Theta) -sind(Theta); sind(Theta) cosd(Theta)];
Sigma = (Sigma*Scale)*Rotate

x1 = -100:1:100; x2 = -100:1:100;
[X1,X2] = meshgrid(x1,x2);
F = mvnpdf([X1(:) X2(:)],mu,Sigma);
F = reshape(F,length(x2),length(x1));

imshow(F*255)

Unfortunately, when I change Theta to a value other than 0, it says SIGMA must be a square, symmetric, positive definite matrix. Can I know what's going on

Answer 1

If you consult the article on Wikipedia about the general elliptical version of the Gaussian 2D PDF, it doesn't look like you're rotating it properly. In general, the equation is:

^{Source: Wikipedia}

where:

Usually, A = 1 and we'll adopt that here. The angle theta will rotate the PDF counter-clockwise, and so we can use this raw form of the equation over mvnpdf. Using your definitions and constants, it would become this:

x1 = -100:1:100; x2 = -100:1:100;
[X1,X2] = meshgrid(X1, X2);
sigma1 = 1;
sigma2 = 1;
scale1 = 3;
scale2 = 1;
sigma1 = scale1*sigma1;
sigma2 = scale2*sigma2;
Theta = 10;

a = ((cosd(Theta)^2) / (2*sigma1^2)) + ((sind(Theta)^2) / (2*sigma2^2));
b = -((sind(2*Theta)) / (4*sigma1^2)) + ((sind(2*Theta)) / (4*sigma2^2));
c = ((sind(Theta)^2) / (2*sigma1^2)) + ((cosd(Theta)^2) / (2*sigma2^2));

mu = [0 0];
A = 1;
F = A*exp(-(a*(X1 - mu(1)).^2 + 2*b*(X1 - mu(1)).*(X2 - mu(2)) + c*(X2 - mu(2)).^2));
imshow(F*255);

Answer 2

I see that you are trying to rotate the variances in 2D... You have to hit your sigma matrix with your transform on both sides of it since the "sigma" matrix is essentially a tensor.

mu = [0 0];

Sigma = [1 0; 0 1];
Scale = [10 0;0 1];
Theta = pi/3;
m = makehgtform('zrotate',Theta);
m = m(1:2,1:2);
Sigma = m*(Sigma*Scale)*m.';

x1 = -100:1:100; 
x2 = -100:1:100;
[X1,X2] = meshgrid(x1,x2);
F = mvnpdf([X1(:) X2(:)],mu,Sigma);
F = reshape(F,length(x2),length(x1));

imshow(F*255)