|
7551
|
1 ## Copyright (C) 2008 Bill Denney |
|
|
2 ## |
|
|
3 ## This file is part of Octave. |
|
|
4 ## |
|
|
5 ## Octave is free software; you can redistribute it and/or modify it |
|
|
6 ## under the terms of the GNU General Public License as published by |
|
|
7 ## the Free Software Foundation; either version 3 of the License, or (at |
|
|
8 ## your option) any later version. |
|
|
9 ## |
|
|
10 ## Octave is distributed in the hope that it will be useful, but |
|
|
11 ## WITHOUT ANY WARRANTY; without even the implied warranty of |
|
|
12 ## MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU |
|
|
13 ## General Public License for more details. |
|
|
14 ## |
|
|
15 ## You should have received a copy of the GNU General Public License |
|
|
16 ## along with Octave; see the file COPYING. If not, see |
|
|
17 ## <http://www.gnu.org/licenses/>. |
|
|
18 |
|
|
19 ## -*- texinfo -*- |
|
|
20 ## @deftypefn {Function File} {@var{area} =} rectint (@var{a}, @var{b}) |
|
|
21 ## |
|
|
22 ## Compute the area of intersection of rectangles in @var{a} and |
|
|
23 ## rectangles in @var{b}. Rectangles are defined as [x y width height] |
|
|
24 ## where x and y are the minimum values of the two orthogonal |
|
|
25 ## dimensions. |
|
|
26 ## |
|
|
27 ## If @var{a} or @var{b} are matrices, then the output, @var{area}, is a |
|
8494
|
28 ## matrix where the i-th row corresponds to the i-th row of a and the j-th |
|
|
29 ## column corresponds to the j-th row of b. |
|
7551
|
30 ## |
|
|
31 ## @seealso{polyarea} |
|
|
32 ## @end deftypefn |
|
|
33 |
|
|
34 ## Author: Bill Denney <bill@denney.ws> |
|
|
35 |
|
|
36 function area = rectint (a, b) |
|
|
37 |
|
|
38 if (nargin != 2) |
|
|
39 print_usage (); |
|
|
40 elseif (ndims (a) != 2 || ndims (b) != 2) |
|
|
41 error ("rectint: expecting arguments to be 2-d arrays"); |
|
|
42 elseif (columns (a) != 4) |
|
|
43 error ("rectint: a must have 4 columns"); |
|
|
44 elseif (columns (b) != 4) |
|
|
45 error ("rectint: b must have 4 columns"); |
|
7655
|
46 elseif any ([a(:,3:4);b(:,3:4)](:) < 0) |
|
8664
|
47 error ("rectint: all widths and heights must be > 0"); |
|
7655
|
48 endif |
|
|
49 |
|
|
50 ## This runs faster if the number of rows of a is greater than the |
|
|
51 ## number of rows of b. Swap them and transpose to make it run |
|
|
52 ## faster. |
|
|
53 swapinputs = false (); |
|
|
54 if (rows (a) > rows (b)) |
|
|
55 tmp = a; |
|
|
56 a = b; |
|
|
57 b = tmp; |
|
|
58 swapinputs = true (); |
|
7551
|
59 endif |
|
|
60 |
|
|
61 area = zeros (rows (a), rows (b)); |
|
7655
|
62 r1 = [a(:,1:2) a(:,1:2)+a(:,3:4)]; |
|
|
63 r2 = [b(:,1:2) b(:,1:2)+b(:,3:4)]; |
|
|
64 for i = 1:columns (area) |
|
|
65 ## Find the location of each point relative to the other points. |
|
|
66 r1x1small = r1(:,1) < r2(i,1); |
|
|
67 r1x1large = r1(:,1) > r2(i,3); |
|
|
68 r1x1mid = (r1(:,1) >= r2(i,1)) & (r1(:,1) <= r2(i,3)); |
|
|
69 r1x2small = r1(:,3) < r2(i,1); |
|
|
70 r1x2large = r1(:,3) > r2(i,3); |
|
|
71 r1x2mid = (r1(:,3) >= r2(i,1)) & (r1(:,3) <= r2(i,3)); |
|
|
72 |
|
|
73 r1y1small = r1(:,2) < r2(i,2); |
|
|
74 r1y1large = r1(:,2) > r2(i,4); |
|
|
75 r1y1mid = (r1(:,2) >= r2(i,2)) & (r1(:,2) <= r2(i,4)); |
|
|
76 r1y2small = r1(:,4) < r2(i,2); |
|
|
77 r1y2large = r1(:,4) > r2(i,4); |
|
|
78 r1y2mid = (r1(:,4) >= r2(i,2)) & (r1(:,4) <= r2(i,4)); |
|
|
79 |
|
|
80 ## determine the width of the rectangle |
|
|
81 ## r1 completely encloses r2 |
|
|
82 area(r1x1small & r1x2large,i) = r2(i,3) - r2(i,1); |
|
|
83 ## the range goes from r2x min to r1x max |
|
|
84 mask = r1x1small & r1x2mid; |
|
|
85 area(mask,i) = r1(mask,3) - r2(i,1); |
|
|
86 ## the range goes from r1x min to r2x max |
|
|
87 mask = r1x1mid & r1x2large; |
|
|
88 area(mask,i) = r2(i,3) - r1(mask,1); |
|
|
89 ## the range goes from r1x min to r1x max |
|
|
90 mask = r1x1mid & r1x2mid; |
|
|
91 area(mask,i) = r1(mask,3) - r1(mask,1); |
|
|
92 |
|
|
93 ## determine the height of the rectangle |
|
|
94 ## r1 completely encloses r2 |
|
|
95 area(r1y1small & r1y2large,i) .*= r2(i,4) - r2(i,2); |
|
|
96 ## the range goes from r2y min to r1y max |
|
|
97 mask = r1y1small & r1y2mid; |
|
|
98 area(mask,i) .*= r1(mask,4) - r2(i,2); |
|
|
99 ## the range goes from r1y min to r2y max |
|
|
100 mask = r1y1mid & r1y2large; |
|
|
101 area(mask,i) .*= r2(i,4) - r1(mask,2); |
|
|
102 ## the range goes from r1x min to r1x max |
|
|
103 mask = r1y1mid & r1y2mid; |
|
|
104 area(mask,i) .*= r1(mask,4) - r1(mask,2); |
|
|
105 |
|
7551
|
106 endfor |
|
|
107 |
|
7655
|
108 if swapinputs |
|
|
109 area = area'; |
|
7551
|
110 endif |
|
|
111 |
|
|
112 endfunction |
|
|
113 |
|
|
114 ## Tests |
|
7655
|
115 ## Exactly overlapping |
|
|
116 %!assert(rectint([0 0 1 1], [0 0 1 1]), 1) |
|
|
117 ## rect2 completely enclosed by rect1 |
|
|
118 %!assert(rectint([-1 -1 3 3], [0 0 1 1]), 1) |
|
|
119 ## rect1 completely enclosed by rect2 |
|
|
120 %!assert(rectint([0 0 1 1], [-1 -1 3 3]), 1) |
|
|
121 ## rect1 right and top in rect2 |
|
|
122 %!assert(rectint([-1 -1 1.5 1.5], [0 0 1 1]), 0.25) |
|
|
123 ## rect2 right and top in rect1 |
|
|
124 %!assert(rectint([0 0 1 1], [-1 -1 1.5 1.5]), 0.25) |
|
|
125 ## no overlap - shared corner |
|
7551
|
126 %!assert(rectint([0 0 1 1], [1 1 2 2]), 0) |
|
7655
|
127 ## no overlap - shared edge |
|
|
128 %!assert(rectint([0 0 1 1], [0 1 2 2]), 0) |
|
|
129 ## Correct orientation of output |
|
|
130 %!assert(rectint([0 0 1 1;0.5 0.5 1 1;-1 -1 2 2], [1 1 2 2]), [0;0.25;0]) |
|
|
131 %!assert(rectint([1 1 2 2], [0 0 1 1;0.5 0.5 1 1;-1 -1 2 2]), [0 0.25 0]) |
