Intersecting graphics

Question

Does the Mathematica graphics system have any concept of intersecting graphics? I've not found much in the documents so far. For example, if I want to show the intersection of two shapes:

Graphics[{Rectangle[], Disk[{0.2, 0}, .5]}]

I know I can use Opacity:

Graphics[{Opacity[0.8], Red, Rectangle[], Green, Disk[{0.2, 0}, .5]}]

But is there a way of specifying the colours of intersecting areas directly? It doesn't seem to be possible to 'address' the intersecting shapes any other way.

In the same vein, is is possible to 'extract' the graphical intersection of arbitrary shapes, without returning to the original geometry and calculating it? Could you obtain this type of entity easily given the above specification (these are just examples...!):

I think it might be easier with raster images, but am interested for now in vector graphics.

Answer 1

How about RegionPlot?

RegionPlot[
  {
   (x - 0.2)^2 + y^2 < 0.5 && 0 < x < 1 && 0 < y < 1,
   (x - 0.2)^2 + y^2 < 0.5 && ! (0 < x < 1 && 0 < y < 1),
   ! ((x - 0.2)^2 + y^2 < 0.5) && 0 < x < 1 && 0 < y < 1
  }, 
   {x, -1, 1.5}, {y, -1, 1.5}, 
   PlotStyle -> {Red, Yellow, Blue}
]

EDIT in response to Szabolcs's comment:

PointInPoly[{x_, y_}, poly_List] := 
 Module[{i, j, c = False, npol = Length[poly]}, 
  For[i = 1; j = npol, i <= npol, j = i++, 
   If[((((poly[[i, 2]] <= y) && (y < 
             poly[[j, 2]])) || ((poly[[j, 2]] <= y) && (y < 
             poly[[i, 2]]))) && (x < (poly[[j, 1]] - 
             poly[[i, 1]])*(y - poly[[i, 2]])/(poly[[j, 2]] - 
              poly[[i, 2]]) + poly[[i, 1]])), c = ¬ c];];
  c]

(from an answer I gave in MathGroup)

RegionPlot[{
   PointInPoly[{x, y}, {{1, 3}, {3, 4}, {4, 7}, {5, -1}, {3, -3}}] && 
   PointInPoly[{x, y}, {{2, 2}, {3, 3}, {4, 2}, {0, 0}}], 
   PointInPoly[{x, y}, {{1, 3}, {3, 4}, {4, 7}, {5, -1}, {3, -3}}] &&
   ¬ PointInPoly[{x, y}, {{2, 2}, {3, 3}, {4, 2}, {0, 0}}],
   ¬ PointInPoly[{x, y}, {{1, 3}, {3, 4}, {4, 7}, {5, -1}, {3, -3}}] &&
   PointInPoly[{x, y}, {{2, 2}, {3, 3}, {4, 2}, {0, 0}}]}, 
  {x, 0, 6}, {y, -4, 8}, 
  PlotPoints -> 100, PlotStyle -> {Red, Yellow, Blue}
]

Answer 2

I'm coming to the party a bit late, but here's my approach. It should work for any two polygons, including non-convex and self-intersecting ones.

winding[poly_, pt_] := 
 Round[(Total@
      Mod[(# - RotateRight[#]) &@(ArcTan @@ (pt - #) & /@ poly), 
       2 Pi, -Pi]/2/Pi)]
cross[e1_, e2_] /; (N[Det[{Subtract @@ e1, Subtract @@ e2}]] === 0.) =
   None;
cross[e1_, e2_] := Module[{params},
   params = ((e2[[2]] - 
        e1[[2]]).Inverse[{Subtract @@ e1, -(Subtract @@ e2)}]);
   If[And @@ Thread[0 <= params <= 1], 
    Subtract @@ e1 params[[1]] + e1[[2]],
    None]];

intersection[poly1_, poly2_, p : {in1_, in2_} : {1, 1}] := 
  Module[{edges1, edges2, intersections,
    inter1, inter2, newedges1, newedges2, midpoints1, midpoints2},
   edges1 = Partition[Range[Length[poly1]], 2, 1, {1, 1}];
   edges2 = Partition[Range[Length[poly2]], 2, 1, {1, 1}];

   intersections = Table[cross[poly1[[e1]], poly2[[e2]]],
     {e1, edges1}, {e2, edges2}];
   inter1 = Flatten[Table[
      SortBy[
       Prepend[DeleteCases[intersections[[i]], None], poly1[[i]]], 
       Norm[# - poly1[[i]]] &], {i, Length[edges1]}], 1];
   inter2 = 
    Flatten[Table[
      SortBy[Prepend[DeleteCases[intersections[[All, i]], None], 
        poly2[[i]]], Norm[# - poly2[[i]]] &], {i, Length[edges2]}], 1];

   newedges1 = Partition[inter1, 2, 1, {1, 1}];
   newedges2 = Partition[inter2, 2, 1, {1, 1}];

   midpoints1 = Mean /@ newedges1;
   midpoints2 = Mean /@ newedges2;
   Flatten[{Pick[newedges1, Abs[winding[poly2, #]] & /@ midpoints1, 
       in1],
      Pick[newedges2, Abs[winding[poly1, #]] & /@ midpoints2, in2]}, 
     1] //.
    {{a___, {b__, c_List}, d___, {c_, e__}, 
       f___} :> {a, {b, c, e}, d, f},
     {a___, {b__, c_List}, d___, {e__, c_}, f___} :> {a, 
       Join[{b, c}, Reverse[{e}]], d, f},
     {a___, {c_List, b__}, d___, {c_, e__}, f___} :> {a, 
       Join[Reverse[{e}], {c, b}], d, f},
     {a___, {c_List, b__}, d___, {e__, c_}, f___} :> {a, {e, c, b}, d,
        f}
     }
   ];

Some notes

winding and cross are two helper functions. winding calculates the winding number of a point pt with respect to a polygon poly given as a list of vertex coordinates. A point lies inside a polygon if and only if the winding number is non-zero.

The function cross calculates the intersection point of two line segments, or returns None if they don't intersect.

intersection is the main function which calculates the intersecting polygon of two polygons poly1 and poly2. It works by calculating the intersection points between the two polygons and adding these to the vertex lists of poly1 and poly2. Then each of the edges of the new polygons lie either completely inside or outside of the other polygon.

The intersection of the two polygons $\text{poly1} \cap \text{poly2}$ is then the union of edges of poly1 that lie inside poly2 and vice versa. Similarly one can also calculate the complement of the two polygons, $\text{poly1} \backslash \text{poly2}$ and $\text{poly1} \backslash \text{poly2}$, and the union $\text{poly1} \cup \text{poly2}$. These four options can be set by in1 and in2.

Example

Manipulate[DynamicModule[{ips11, ips10, ips01},
  pts = PadRight[pts, 2 n, RandomReal[{-1, 1}, {2 n, 2}]];
  ips11 = intersection[pts[[ ;; n]], pts[[n + 1 ;;]], {1, 1}];
  ips10 = intersection[pts[[ ;; n]], pts[[n + 1 ;;]], {1, 0}];
  ips01 = intersection[pts[[ ;; n]], pts[[n + 1 ;;]], {0, 1}];
  Graphics[{
    {Yellow, Polygon[ips10]},
    {Blue, Polygon[ips01]},
    {Red, Polygon[ips11]},
    {FaceForm[], EdgeForm[Black], Polygon[pts[[ ;; n]]]}, {FaceForm[], 
     EdgeForm[Black], Polygon[pts[[n + 1 ;;]]]}}, 
   PlotRange -> {{-1, 1}, {-1, 1}}]], 
  {{pts, {}}, Locator}, {{n, 5}, None}]

Answer 3

The (undocumented!) function Graphics`PolygonUtils`PolygonIntersection[] (Graphics`Mesh`PolygonIntersection[] in older versions) seems up to the task. Using Sjoerd's example:

polys = {Polygon[{{1, 3}, {3, 4}, {4, 7}, {5, -1}, {3, -3}}],
         Polygon[{{2, 2}, {3, 3}, {4, 2}, {0, 0}}]};

Graphics[Append[{Gray, polys}, {Blue, Graphics`PolygonUtils`PolygonIntersection[polys]}]]

Disk[] objects are not covered by this method, but it is not too difficult to make a Polygon[] that approximates a Disk[]...

Answer 4

I am not aware of any built-in functionality (I might easily be wrong), but there's an example at MathWorld for calculating intersections of convex polygons. You'd need to approximate the circle with a polygon.

Get the notebook from that page: there's an intersection calculation inside that uses the IMTEK Mathematica supplement.

Example:

<< Imtek`Polygon`

disk = Disk[{0.2, 0}, 0.5];
rec = Polygon[{{0, 0}, {1, 0}, {1, 1}, {0, 1}}];

Graphics[{Green, rec, Red, disk, Blue, 
  Polygon@imsConvexIntersect[{imsPolygonizeCircle[Circle @@ disk, 50],
      First[rec]}]}]

Answer 5

Another option is to use image processing features such as ImageCompose:

{
 g1 = Graphics[Rectangle[], PlotRange -> 1],
 g2 = Graphics[Disk[{0.2, 0}, .5], PlotRange -> 1],
 ImageCompose[g1, g2, Center, Center, {1, 1, 0}]
 }

The output of the above looks like this:

(Note that in this case your Graphics get rasterized and the result is an Image.)

Answer 6

In version 10, the new geometric computation functionality supports this. It operates with region objects. Many graphics primitives, including Disk and Rectangle can be used as regions.

Boolean operations include RegionUnion, RegionIntersect, RegionDifference, RegionSymmetricDifference and BooleanRegion.

Example:

RegionPlot[RegionIntersection[Rectangle[], Disk[{0.2, 0}, .5]], AspectRatio->Automatic]

In many cases these Boolean operators do not evaluate. RegionIntersection[Rectangle[], Disk[{0.2, 0}, .5]] can be used as a region, but it does not evaluate to another expression. However, it is possible to get an approximation to the result, free of any explicit Boolean expressions, using DiscretizeRegion. It is then possible to extract coordinates (e.g. the boundary) from this discretized mesh.

Answer 7

You may be able to do this using FilledCurve in version 8. I see examples of subtraction and exclusion but not intersection.