Code covered by the BSD License

### Highlights fromgeom2d

3.86667

3.9 | 18 ratings Rate this file 281 Downloads (last 30 days) File Size: 221 KB File ID: #7844

# geom2d

13 Jun 2005 (Updated 29 Feb 2012)

Geometry library for matlab. Performs geometric computations on points, lines, circles, polygons...

File Information
Description

Geometry library to handle and visualize geometric primitives such as points, lines, circles and ellipses, polylines and polygons...

The goal is to provide a low-level library for manipulating geometrical primitives, making easier the development of more complex geometric algorithms.

The library proposes functions to:

- create various shapes (points, circles, lines, ellipses, polylines, polygons) using an intuitive syntax. Ex: createCircle(p1, p2, p3) to create a circle through 3 points.

- derive new shapes: intersection between 2 lines, between a line and a circle, parallel and perpendicular lines
Functions to compute intersections

- work on polylines and polygons: compute centroid and area, expand, clip with half-plane...

- measure distances (between points, a point and a line, a point and a group of points), angle (of a line, between 3 points), or test geometry (point on a line, on a circle).

- manipulate planar transformation. Ex: P2 = transformPoint(P1, createRotation(CENTER, THETA));

- draw shapes easily. Ex: drawCircle([50 50], 25), drawLine([X0 Y0 DX DY]). Some clipping is performed for infinite shapes such as lines or rays.

Additional help is provided in geom/Contents.m file, as well as summary files like 'points2d.m' or 'lines2d.m'.

Note: the project has merged with the geom3d library (FeX 24484), and is now hosted on sourceforge: http://matgeom.sourceforge.net/

Acknowledgements

This file inspired Image Ellipsoid 3 D.

MATLAB release MATLAB 7.9 (R2009b)
Other requirements Some functions (for shape fitting) require the optimization toolbox.
Tags for This File
Everyone's Tags
Tags I've Applied
08 Mar 2013

Really nice works... Thanks

08 Mar 2012

Hi sumitha,
basic format for representing lines is [x0 y0 dx dy], where (w0,y0) refers to origin point, and (dx,dy) to direction vector. You can create a line orthogonal to direction THETA with:
line1 = [x0 Y0 -sin(THETA) cos(THETA)];
drawLine(line1);
Regards.

07 Mar 2012

Hi,can u please suggest me which one of these i can use to draw a line orthogonal to the given orientation(in radian)or else can u plz give me an idea how i can proceed with the same.

01 Mar 2012

Hi Benjamin, the clipping of polygon is not an easy problem. There are efficient algorithms, and their inclusion is planned, but I still haven't found time to write them down. It seems the mapping toolbox provides facilities for this, maybe you can check it.
Regards, David

28 Feb 2012

Hi David, thanks for the great work. I have been using your clipPolygon successfully for clipping convex polygons, but now I also need it for non-convex polygons... are you planning to include this in a future release?

18 Feb 2012

Very nice and essential coding thanks for sharing.

19 Nov 2011

06 Jun 2011

@Frank: yes, there was a stupid bug in this function. I have submitted an enhanced version that fixes it.
Regards, David

02 Jun 2011

Hi, looks like a nice little toolbox!
I just tried your function "intersectLineCircle.m". Unfortunately, it seems that the code does not give the right answers for tangential lines (Returns 'nan') and lines which are not intersecting at all (gives a wrong intersection point).

Here is the data I used:
c1 = [0 0 0.5]
l1 = [0.5 0 0 1] %tangential
l1 = [-1 0 0 1] %non touching
pts = intersectLineCircle(l1, c1)

Cheers,
Frank

08 Apr 2011

geom2d is a first class complement to the geometry tools provided.
David: You should add a 2D alpha shape routine. I am testing this one at the moment, but leaves lots of room for improvement. http://www.mathworks.com/matlabcentral/fileexchange/6760-ashape-a-pedestrian-alpha-shape-extractor

16 Dec 2010

Hi Chris,
intersection of circle with line/circle is not implemented for the moment, juste because I did not need it yet... But I will try to implement it in a future release.
Regards.

11 Dec 2010

Hi, this is very useful for me? I just wonder why didn't i find the function for obtaining intersections between a line and a circle? Thanks.

01 Dec 2010

Well done, David!

14 Sep 2010

OK, I think this should work, circumventing eps. Again, looking at center points instead of end-points. This time the center points of segments. It still needed some fail-safe code when multiple edge segments show up inside box due to finite precission. The code should also maintain 'direction' of the edges, as your previous code.
Here the conplete code for clipEdge:

% process data input
if size(box, 1)==1
box = box([1 3 2 4]);
end

% get limits of window
xmin = box(1);
ymin = box(2);
xmax = box(3);
ymax = box(4);

% convert window limits into lines
lineX0 = [xmin ymin xmax-xmin 0];
lineX1 = [xmin ymax xmax-xmin 0];
lineY0 = [xmin ymin 0 ymax-ymin];
lineY1 = [xmax ymin 0 ymax-ymin];

edge2 = zeros(size(edge));

for i=1:size(edge,1)
% current edge
iedge = edge(i, :);

% compute intersection points with each line of bounding window
px0 = intersectLineEdge(lineX0, iedge);
px1 = intersectLineEdge(lineX1, iedge);
py0 = intersectLineEdge(lineY0, iedge);
py1 = intersectLineEdge(lineY1, iedge);

% create array of points
points = [px0; px1; py0; py1; iedge(1:2); iedge(3:4)];

% and remove infinite points (edges parallel to box edges)
points = points(all(isfinite(points),2), :);

% sort points by x then y
points = sortrows(points);

% get center positions between points
centerPoints=(points(2:end,:)+points(1:end-1,:))/2;

% find the centerPoints (if any) inside window
ins = find( centerPoints(:,1)>=xmin & centerPoints(:,2)>=ymin & ...
centerPoints(:,1)<=xmax & centerPoints(:,2)<=ymax);

% if multiple segments inside box, which can happen due to finite resolution, only take the longest segment
if length(ins)>1
dv=points(ins+1,:)-points(ins,:); % delta vectors of the segments
len=hypot(dv(:,1), dv(:,2)); % lengths of segments
[a,I]=max(len); % find index of longest segment
ins=ins(I);
end

if length(ins)==1 % if one of the center points is inside box, then the according edge segment is indide box
if edge(i,1)>edge(i,3) || (edge(i,1)==edge(i,3) && edge(i,2)>edge(i,4)) % restore same direction of edge
edge2(i, :) = [points(ins+1,:) points(ins,:)];
else
edge2(i, :) = [points(ins,:) points(ins+1,:)];
end
end
end

14 Sep 2010

Argh, never mind. My solution has div by 0 problem. The whole eps stuff is really a pain in the ...

14 Sep 2010

turns out my suggestion above introduces a problem with negative numbers. The code below should work correctly:

% ins = find( points(:,1)-xmin>=-eps & points(:,2)-ymin>=-eps & ...
% points(:,1)-xmax<=eps & points(:,2)-ymax<=eps);
ins = find( (points(:,1)/xmin-1)*sign(xmin)>=-eps & (points(:,2)/ymin-1)*sign(ymin)>=-eps & ...
(points(:,1)/xmax-1)*sign(xmax)<=eps & (points(:,2)/ymax-1)*sign(ymax)<=eps);

13 Sep 2010

Also in clipEdge() you should change the limit check. See code below. Change commented out to suggested.

cheers,
Reto

% ins = find( points(:,1)-xmin>=-eps & points(:,2)-ymin>=-eps & ...
% points(:,1)-xmax<=eps & points(:,2)-ymax<=eps);
ins = find( points(:,1)/xmin-1>=-eps & points(:,2)/ymin-1>=-eps & ...
points(:,1)/xmax-1<=eps & points(:,2)/ymax-1<=eps);

09 Sep 2010

Bonjour,
again some trouble with large numbers and limit checking. This time in isPointOnEdge, which is used for intersectLineEdge, which in turn is used in clipEdges.
isPointOnEdge will not work correctly for coordinates with large numbers, as you apply a fixed tolerance. However, the numerical tolerance of matlab is relative to the number itself.
I normalized the point and edge x and y elements to the point x and y respectively. Like that you will deal with numbers in the range of 1 and your tolerance check will work.
See the code below. Note that I also dumped the Np==Ne check, as it didn't do anything.
Regards, Reto

% extract computation tolerance
tol = 1e-14;
if ~isempty(varargin)
tol = varargin{1};
end

% number of edges and of points
Np = size(point, 1);
Ne = size(edge, 1);

x0 = repmat(edge(:,1)', Np, 1);
y0 = repmat(edge(:,2)', Np, 1);
dx = repmat(edge(:,3)', Np, 1)-x0;
dy = repmat(edge(:,4)', Np, 1)-y0;
xp = repmat(point(:,1), 1, Ne);
yp = repmat(point(:,2), 1, Ne);

% normalize x, y to xp, yp
x0 = x0/xp;
y0 = y0/yp;
dx = dx/xp;
dy = dy/yp;
xp = 1;
yp = 1;

% test if point is located on supporting line
b1 = abs((xp-x0).*dy - (yp-y0).*dx)./hypot(dx, dy)<tol;

% check if point is located between edge bounds
% use different tests depending on line angle
ind = abs(dx)>abs(dy);
t = zeros(max(Np, Ne), 1);
t(ind) = (xp(ind)-x0(ind))./dx(ind);
t(~ind) = (yp(~ind)-y0(~ind))./dy(~ind);
b = t>-tol & t-1<tol & b1;

05 Aug 2010

Hi,
this is a good idea, much more convenient than dealing with eps. I will implement it soon.

04 Aug 2010

Hi David,
your new limit check assumes that the direction vector of the line has a similar scale as the bounding box, which is not necessarely true. The direction vector often has a magnitude of 1, while the bounding box can be numerically very large.

Here is an idea and the thinking behind it:
You are checking the endpoints of the edge to be within the bounding box. Naturally, the end points are always on the bounding box, so you run into numerical precission issues. So why not check if the center of the edge is within the box? This way you are away from the limit lines with the exception of horizontal and vertical lines. In those two cases you should not have any numerical issues, as they are at right angles.

So, I would re-write the limit check as follows and drop eps completely:

% check if edge center is within bounds
edgeCenterX=sum(edge(i,[1,3]))/2;
xMinOk = xmin <= edgeCenterX;
xMaxOk = edgeCenterX <= xmax;
edgeCenterY=sum(edge(i,[2,4]))/2;
yMinOk = ymin <= edgeCenterY;
yMaxOk = edgeCenterY <= ymax;

It works for me.

04 Aug 2010

Hi Reto,

thank you for your numerous and constructive remarks. I have published a new version that take into account your corrections. I will try to propagate modification to other drawing functions.

For limit check, I have used a different strategy, but it seems to work with large bouding boxes as well. Tell me if there are still problems.

03 Aug 2010

Also consider re-writing adding the precision to the limit check as follows:

xOk = xmin*(1-eps)-min(edge(:,[1,3]), [], 2)<0 & max(edge(:,[1,3]), [], 2)-xmax*(1+eps)<0;
yOk = ymin*(1-eps)-min(edge(:,[2,4]), [], 2)<0 & max(edge(:,[2,4]), [], 2)-ymax*(1+eps)<0;

This will take into consideration the scale of the limiting box. With the original code you run into precission issues when the box has large numbers (I'm working with numbers in the 1e9 range). Like this the precision is relative rather than absolute.

03 Aug 2010

Yet an other bug. This library is a good start, but definitely not ready fro prime time.
in clipLine() there is incorrect indexing. Fix is below. Consider using this fixed function in drawLine() (see bug in drawLine() above).

% check bounds of result edge in each direction

% incorrect indices, indexing a edge, not a box
% xOk = xmin-min(edge(:,1:2), [], 2)<eps & max(edge(:,1:2), [],2)-xmax<eps;
% yOk = ymin-min(edge(:,3:4), [], 2)<eps & max(edge(:,3:4), [],2)-ymax<eps;

% this is correct indexing
xOk = xmin-min(edge(:,[1,3]), [], 2)<eps & max(edge(:,[1,3]), [], 2)-xmax<eps;
yOk = ymin-min(edge(:,[2,4]), [], 2)<eps & max(edge(:,[2,4]), [], 2)-ymax<eps;

02 Aug 2010

oops, now my correction of the correction replaced the original correction. You can tell I'm a novice to this system...
So, here again the complete suggested correction to drawLine().
The original code will draw horizontal lines even though they are outside the axis. It will also draw lines in the quadrant to the right and below the axis.

% sort points along the x coordinate, and draw a line between
% the two in the middle
points = sortrows([px1 ; px2 ; py1 ; py2], 1);

if isfinite(points(3,1)) % if not horizontal or vertical
points=points(2:3,:);
else % else horizontal or vertical
points=points(1:2,:);
end

if points(1,1)>=xmin && points(2,1)<=xmax &&...
points(1,2)>=ymin && points(2,2)<=ymax

h(i)=line(points(1:2, 1), points(1:2, 2));

if ~isempty(varargin)
set(h(i), varargin{:});
end

else
h(i)=-1;
end

02 Aug 2010

Sorry, correction to the correction:

if points(1,1)>=xmin && points(2,1)<=xmax &&...
points(1,2)>=ymin && points(2,2)<=ymax

02 Aug 2010

There is a bug in intersectEdges(). Change lines 143 and 144 from

x0(out) = NaN;
x0(out) = NaN;

to

x0(i(out)) = NaN;
y0(i(out)) = NaN;

otherwise indexing is incorrect and you get false results.

20 Jul 2010

@Matthew:
I have submitted a new version, that fixes the bug in the 'centroid' function, as well as many others.
regards

13 Jul 2010

Minor (but fatal) bug found in centroid.m for calculating mass-weighted centroid of 2d object.

Change line 69 of the file from:

center = mean(pts.*repmat(mass, 1, 3));

to

center = mean(pts.*repmat(mass, 1, 2));

Or it will error with a subscript assignment dimension mismatch...pts is a 2-column matrix, not a 3-column matrix.

03 Jun 2010
19 May 2010

Hi Aleksandar,
at the moment there is no additional argument for drawRect function.
To specify the figure to draw in, the best is to use figure(...) before calling drawRect, eventually using hold on to avoid clearing previous image.
To change rectangle color or width, you can use the handle returned by the function:
hr = drawRect([x0 y0 w h]);
set(hr, 'color', 'r', 'linewidth', 4);

17 May 2010

I need help with function drawRect.
Can this function have more arguments? (Like: image on which to put rectangle,rectangle color ...)

I have this part of code
im = drawRect(...
picture,...
round((vj(pos, 2) - 1) * sF + 1),...
round((vj(pos, 2) + velicinaLica(1) - 1) * sF + 1),...
round((vj(pos, 3) - 1) * sF + 1),...
round((vj(pos, 3) + velicinaLica(2) - 1) * sF + 1),...
[1 0 0]...
);

24 Jul 2009

oops, sorry for badly corrected bug...
A new version has been submitted, with suggested corrections. There are also some new functions for polygons and polylines (see for example "help polygons2d" or "help polylines2d")

23 Jul 2009

I tested some of the other functions. Quite useful indeed! Thanks.

22 Jul 2009

The file intersectLinePolygon.m still contains 'pi' in place of 'intersects'. The HISTORY says it has been replaced, but it has NOT been replaced yet. I downloaded the zip file on July 22, 2009.

Also, here is a clearer (and of course equivalent way) to do the update:

intersects = [];
for i=1:length(edges)
p0 = intersectLines(line, createLine(edges(i,1:2), edges(i,3:4)));
if isPointOnEdge(p0, edges(i,:))
intersects = [intersects; p0];
end
end

17 Jul 2009

This library and its companion, geom3d, are shining jewels - indispensable for doing non-trivial calculations on objects in 2d and 3d space.

By the way, leave all the files in the geom2d directory and add the directory to your path. Then you can type "help geom2d" or "doc geom2d" and get properly linked help text in the command window or the help window. The Contents.m file is the proper "MATLAB way" of providing documentation.

18 Jun 2009

Hi,

I've submitted a new version, which fixes the bug in intersectEdges, and adds new functions for manipulating polylines (intersections, distance to point...). Please do not hesitate to report bugs (direct email preferred).

The 'pi' variable (cf. alexia's comment) was not supposed to cause problem, as it is located inside of a function. The goal of the library is mainly to hide computational details into functions with explicit names, without having to read the inner code. However, the variable 'pi' was renamed in a more appropriate (and hopefully explicit) name.

regards, David

11 Jun 2009

I've tried the function intersectEdges but found that some cases are mishanlded etc.

23 Nov 2008

Hi, what a bad name for a variable/function "pi", I mean:

"pi=intersectLinePolygon(line, points)"

I had a headache trying to discover what was wrong with this file (because in theory it should work), it was part of a very large project I was working in.

Finally I discovered the "pi" issue, "pi" is present in the code, but I also had to work with pi=3.1415..., and this cr*p of matlab got confused. Oh Mathematica !, I have such good memories...

Certainly is not so clever call a variable/function "pi" in a geometry script. Besides, as peter said is very badly commented (found errors), so don't trust to much with this one.

23 Nov 2008
06 Jun 2007

very helpful in implementing more complicated geometry algorithms

08 Mar 2007

drawLabels doesn't actually draw any labels...

07 Nov 2005

Tried circleaspolygon and rectaspolygon but both were not self contained, badly explained. The package seems to be promising but far from plug and play.

14 Jun 2005

Very handy functions for geometric analysis

12 Apr 2006

update doc, add some functions, fix some bugs

06 Jun 2006

update description

05 Nov 2010

bug fixes (isPointOnEdge, intersectPolylines), remove obsolete licence info

05 Nov 2010

fix bugs (intersectLines, isPointOnEdge), remove obsolete licence considerations

30 Mar 2011

split library in several packages, add several functions (circle and circle-line intersections, computation of inertia ellipse)

06 Jun 2011

fixed bugs in intersectLineCircle, added minimumCaliperDiameter and averagePointSet of a set of points. Now uses degrees for representing shapes

29 Feb 2012

lot of code cleanup, new functions for working with polygons, clean up drawing functions