CODEWITHSUNDEEP
matlabmatrix
Pedro77
Pedro77

Reputation: 5294

Best way to get the bounding rectangle of a set of 3D points on a plane in Matlab

I need to get the four edges of the bounding rectangle of a set of 3D points stored as a 3xN matrix (tt). N >=4. The points lies on a plane.

Code sample:

% Simulate some points
deltaXY = 20;
[xx,yy] = meshgrid(-100:deltaXY:100,-100:deltaXY:100);
XYZ = [xx(:)'; yy(:)'; zeros(1,numel(xx))];

% Add some imperfection to data removing the top rigth point
maxXids = find(XYZ(1,:) == max(XYZ(1,:)));
maxYids = find(XYZ(2,:) == max(XYZ(2,:)));
id = intersect(maxXids,maxYids);
XYZ = removerows(XYZ',id)';

% Lets rotate a bit
XYZ = roty(5)*rotx(7)*rotz(0)*XYZ;

% Plot points
figure;
grid on;
rotate3d on;
axis vis3d;
hold on;
plot3(XYZ(1,:),XYZ(2,:),XYZ(3,:),'.r');

% Find bounding rectangle
% ??? :(
%Currently I'm using this code:
tt = XYZ;
%Get the max and min indexes
minX = find(tt(1,:) == min(tt(1,:)));
minY = find(tt(2,:) == min(tt(2,:)));
maxX = find(tt(1,:) == max(tt(1,:)));
maxY = find(tt(2,:) == max(tt(2,:)));
%Intersect to find the common index
id1 = intersect(minX,minY);
id2 = intersect(maxX,minY);
id3 = intersect(maxX,maxY);
id4 = intersect(minX,maxY);
%Get the points
p1 = tt(:,id1(1));
p2 = tt(:,id2(1));
p3 = tt(:,id3(1));
p4 = tt(:,id4(1));

Sample points plot: Sample points plot

The problem is that intersect some times can be null, eg: if the points does not form a rectangle. Resulting this error:

Index exceeds matrix dimensions.

Upvotes: 1

Views: 584

Answers (2)

Pedro77
Pedro77

Reputation: 5294

Ok found a solution:

% Find bounding rectangle
tt = XYZ;

%Get the max and min indexes
minXids = find(tt(1,:) == min(tt(1,:)));
minYids = find(tt(2,:) == min(tt(2,:)));
maxXids = find(tt(1,:) == max(tt(1,:)));
maxYids = find(tt(2,:) == max(tt(2,:)));

%Intersect to find the common index
id1 = intersect(minXids,minYids);
id2 = intersect(maxXids,minYids);
id3 = intersect(maxXids,maxYids);
id4 = intersect(minXids,maxYids);

%Get the points
% Find affine plane on points
[np,~,pp] = affine_fit(tt');
% Converts to cart. eq.
% ax + yb + cz + d = 0
% Find d
a = np(1); b = np(2); c = np(3);
x = pp(1); y = pp(2); z = pp(3);
d = - (a*x + y*b + c*z);

% Get only one value
minX = min(tt(1,minXids)); maxX = max(tt(1,maxXids));
minY = min(tt(2,minYids)); maxY = max(tt(2,maxYids));

if numel(id1) == 0
    x = minX; y = minY;
    % Calculate z at xy.
    z = - (d + a*x + y*b)/c;
    p1 = [x y z]';
else
    p1 = tt(:,id1(1));
end
if numel(id2) == 0
    x = maxX; y = minY;
    z = - (d + a*x + y*b)/c;
    p2 = [x y z]';
else
    p2 = tt(:,id1(1));
end
if numel(id3) == 0
    x = maxX; y = maxY;
    z = - (d + a*x + y*b)/c;
    p3 = [x y z]';
else
    p3 = tt(:,id1(1));
end
if numel(id4) == 0
    x = minX; y = maxY;
    z = - (d + a*x + y*b)/c;
    p4 = [x y z]';
else
    p4 = tt(:,id1(1));
end

ps = [p1 p2 p3 p4];

plot3(ps(1,:),ps(2,:),ps(3,:),'ob');

Upvotes: 0

BillBokeey
BillBokeey

Reputation: 3476

First solution : Use logical indexing to get rid of the find calls 

p1=tt(:,tt(1,:)==min(tt(1,:))&tt(2,:)==min(tt(2,:)));
p2=tt(:,tt(1,:)==max(tt(1,:))&tt(2,:)==min(tt(2,:)));
p3=tt(:,tt(1,:)==max(tt(1,:))&tt(2,:)==max(tt(2,:)));
p4=tt(:,tt(1,:)==min(tt(1,:))&tt(2,:)==max(tt(2,:)));

Second solution : Use convhull to get the corners :

k=convhull(tt(1,:),tt(2,:));
Corners=[tt(:,k(1:end-1))];

Upvotes: 1

Related Questions