检查点是否在多维度椭圆体内

瓦斯克

我在确定某些点是否位于多维椭圆体内部时遇到问题。真正的问题是由潜在的转换(旋转,平移)引起的,我不知道如何将其应用于我的计算中。

在我的工作中,我将找到一些点,需要找出包围椭圆体的最小体积。然后,我将不得不检查其他点集,并确定其中哪些点属于找到的椭球。

我决定从这里使用代码

function [] = Checkup()
points  = [[ 0.53135758, -0.25818091, -0.32382715] 
    [ 0.58368177, -0.3286576,  -0.23854156,] 
    [ 0.18741533,  0.03066228, -0.94294771] 
    [ 0.65685862, -0.09220681, -0.60347573]
    [ 0.63137604, -0.22978685, -0.27479238]
    [ 0.59683195, -0.15111101, -0.40536606]
    [ 0.68646128,  0.0046802,  -0.68407367]
    [ 0.62311759,  0.0101013,  -0.75863324]];
P = points\'; % <- remove the \ symbol here
dimen = length(P(:,1));

% c - vector with ellipse centers
[A, c] = MinVolEllipse(P, 0.01);
[~, Q, V] = svd(A);
radiuses = 1:dimen;

% Calculate radiuses
for i = 1:dimen
    radiuses(1, i) = 1 / sqrt(Q(i,i));
end

% Check if points lie within ellipse and print Ok for every point inside ellipsoid
for i = 1:length(P(1,:)) % length(P(1,:)) is number of points
    value = 0;
    for j = 1:dimen
        %adding ((p_i - c_i) / (r_i))^2 value
        value = value + ( ( (P(j,i) - c(j, 1) )^2 ) / (radiuses(1, j)^2));

    end
    if value <= 1
        disp('Ok')
    end
end

目前,此代码不会打印确定文本(但应将其打印8次)。

根据:“V是给你的椭球的方向旋转矩阵”我认为这是我需要用我的代码工作的最后一个元素。

我的方法

好的,所以我对代码进行了一些更改。现在,我尝试执行以下操作:

假设我有椭球的中心及其半径。

对于每一点:

  1. 更新位置,如下所示:point = point-center->换句话说,将其转换为椭圆体使其居中于点(0,0,... 0)

  2. 像椭圆体平行于所有轴一样旋转:point = invert(V)* point

  3. 使用以下公式检查点是否在椭球内:(point.x / radius.x)^ 2 + ... +(point.i / radius.i)^ 2其中i是维数

  4. 如果方程式给出结果<= 1,则点位于椭圆内。

这种方法好吗?我的意思是-我不知道是否允许反转V矩阵并像将先前的旋转反转一样使用它...

正确的解决方案

看来我提出的方法不错,但还有更好的方法:

function [result] = Ellipse()
% Returns vector consisting of 10 entries
% that represent error ratio for every test set 

result = 0:9;

% Run tests for all point sets
for pointSet = 0:9
    [Ptraining, Ptest] = LoadPointSet(pointSet);
    count = 0;

    % c - vector with ellipse centers
    [A, c] = MinVolEllipse(Ptraining, 0.00001);

    % Check if points lie within ellipse
    for i = 1:length(Ptest(1,:)) % length(P(1,:)) is number of points
        pp = Ptest(:, i) - c;
        value = pp' * A * pp;    
        if value > 1
            count = count + 1;
        end
    end

    % Get the error ratio
    index = pointSet + 1;
    result(index) = count / length(Ptest(1,:));
end

end
大卫·汉门

目前,此代码不会打印确定文本(但应将其打印8次)。

首先,没有理由期望MinVolEllipse函数之类的近似算法会给出包围椭圆的确切最小体积。您正在提供的容差MinVolEllipse,因此您应该预期该函数的结果会有些错误。

更重要的是,您无需在此处使用奇异值分解(但请参阅下文了解详细信息)。椭球表面的方程为(xc)T A(xc)= 1。您需要检查的是每个点的(xc)T A(xc)是否小于1(加上一些公差):

tol_mvee = 0.01;
tol_dist = 0.1;

[A, c] = MinVolEllipse(P, tol_mvee);

% Check if points lie within ellipse and print Ok for every point inside ellipsoid
for i = 1:length(P(1,:)) % length(P(1,:)) is number of points
    x   = P(i,:);
    x_c = x-c;
    d   = dot(x_c, A*x_c);
    if d < 1+tol_dist
        disp('Ok');
    end
end

您的第一点显然是在椭圆体内。所有其他七个点都在包围椭圆体的真实最小体积的边界上或非常接近该边界。该算法将对此产生一些麻烦,并且由于生成的椭圆体明显是椭圆体而不是球形(请参见下面的条件编号),因此也会有一些麻烦。

奇异值分解可能对告诉您A来自矩阵MinVolEllipse是否病态很有用,在此特定问题中就是这种情况。矩阵的条件数超过200,这意味着您最好不要期望得到容差远小于1e-6的结果。

本文收集自互联网,转载请注明来源。

如有侵权,请联系 [email protected] 删除。

编辑于
0

我来说两句

0 条评论
登录 后参与评论

相关文章