:

Triangulation, hatching and texturing of simple polygons

Maple

Hi,

As an amusement,  I decided several months ago to develop some procedures to fill a simple polygon* by hatches or simple textures.

* A simple polygon is a polygon  whose sides either do not intersect or either have a common vertex.

This work started with the very simple observation that if I was capable to hatch or texture a triangle, I will be too to hatch or texture any simple polygon once triangulated.
I also did some work to extend this work to non-simple polygons but there remains some issues to fix (which explains while it is not deliverd here).

The main ideat I used for hatching and texturing is based upon the description of each triangles by a set of 3 inequalities that any interior point must verify.
A hatch of this triangle is thius a segment whose each point is interior.
The closely related idea is used for texturing. Given a simple polygon, periodically replicated to form the texture, the set of points of each replicate that are interior to a given triangle must verify a set of inequalities (the 3 that describe the triangle, plus N if the pattern of the texture is a simple polygon with N sides).

Unfortunately I never finalise this work.
Recently @Christian Wolinski asked a question about texturing that reminded me this "ancient" work of mine.
So I decided to post it as it is, programatically imperfect, lengthy to run, and most of all french-written for a large part.
I guess it is a quite unorthodox way to proceed but some here could be interested by this work to take it back and improve it.

The module named "trianguler" (= triangulate) is a translation into Maple of Frederic Legrand's Python code (full reference given in the worksheet).
I added my own procedure "hachurer" (= hatching) to this module.
The texturing part is not included in this module for it was still in development.

A lot of improvements can be done that I could list, but I prefer not to be too intrusive in your evaluation of this work. So make your own idea about it and do not hesitate to ask me any informations you might need (in particular translation questions).

PS: this work has been done with Maple 2015.2

 > restart:

Reference: http://www.f-legrand.fr/scidoc/docmml/graphie/geometrie/polygone/polygone.html
(in french)
reference herein : M. de Berg, O. Cheong, M. van Kreveld, M. Overmars,
Computational geometry,  (Springer, 2010)

Direct translation of the Frederic Legrand's Python code

Meaning of the different french terms

voisin_sommet  (n, i, di)
let L the list [1, ..., n] where n is the number of vertices
This procedure returns the index of the neighbour (voisin) of the vertex (sommet) i when L is rotated by di

equation_droite  (P0, P1, M)
Let P0 and P1 two vertices and M an arbitrary point.
Let (P0, P1) the vector originated at P0 and ending at P1 (idem for (P0, M)) and  the unitary vector in the Z direction.
This procedure returns (P0, P1) o (P0, M) .

point_dans_triangle  (triangle, M) P1, P2]
This procedure returns "true" if point M is within (strictly) the  triangle "triangle" and "false" if not.

sommet_distance_maximale  (polygone, P0, P1, P2, indices)
Given a polygon (polygone) threes vertices P0, P1 and P2 and a list of indices , this procedure returns
the vertex of the polygon "polygone" which verifies: 1/ this vertex is strictly within
the triangle [P0, P1, P2] and 2/ it is the farthest from side [P1, P2] amid all the vertices that verifies point 1/.
If there is no such vertex the procedure returns NULL.

sommet_gauche (polygone)
This procedure returns the index of the leftmost ("gauche" means "left") vertex in polygon "polygone".
If more than one vertices have the same minimum abscissa value then only the first one is returned.

nouveau_polygone(polygone,i_debut,i_fin)
This procedure creates a new polygon from index i_debut (could be translated by i_first) to i_end (i_last)

trianguler_polygone_recursif(polygone)
This procedure recursively divides a polygon in two parts A and B from its leftmost vertex.
If A (B) is a triangle the list "liste_triangles" (mening "list of triangles") is augmented by A (B);
if not the procedure executes recursively on A and B.

trianguler_polygone(polygone)
This procedure triangulates the polygon "polygon"

hachurer(shapes, hatch_angle, hatch_number, hatch_color)
This procedure generates stes of hatches of different angles, colors and numbers

Limitations:
1/ "polygone" is a simply connected polygon
2/  two different sides S and S', either do not intersect or either have a common vertex

 > trianguler := module() export voisin_sommet, equation_droite, interieur_forme, point_dans_triangle, sommet_distance_maximale,        sommet_gauche, nouveau_polygone, trianguler_polygone_recursif, trianguler_polygone, hachurer: #------------------------------------------------------------------- voisin_sommet := (n, i, di) -> ListTools:-Rotate([\$1..n], di)[i]: #------------------------------------------------------------------- equation_droite := proc(P0, P1, M) (P1[1]-P0[1])*(M[2]-P0[2]) - (P1[2]-P0[2])*(M[1]-P0[1]) end proc: #------------------------------------------------------------------- interieur_forme := proc(forme, M)   local N:   N := numelems(forme);   { seq( equation_droite(forme[n], forme[piecewise(n=N, 1, n+1)], M) >= 0, n=1..N) } end proc: #------------------------------------------------------------------- point_dans_triangle := proc(triangle, M)   `and`(           is( equation_droite(triangle[1], triangle[2], M) > 0 ),           is( equation_droite(triangle[2], triangle[3], M) > 0 ),           is( equation_droite(triangle[3], triangle[1], M) > 0 )        ) end proc: #------------------------------------------------------------------- sommet_distance_maximale := proc(polygone, P0, P1, P2, indices)   local n, distance, j, i, M, d;   n        := numelems(polygone):   distance := 0:   j        := NULL:   for i from 1 to n do     if `not`(member(i, indices)) then       M := polygone[i];       if point_dans_triangle([P0, P1, P2], M) then         d := abs(equation_droite(P1, P2, M)):         if d > distance then           distance := d:           j := i         end if:       end if:     end if:   end do:   return j: end proc: #------------------------------------------------------------------- sommet_gauche := polygone -> sort(polygone, key=(x->x[1]), output=permutation)[1]: #------------------------------------------------------------------- nouveau_polygone := proc(polygone, i_debut, i_fin)   local n, p, i:   n := numelems(polygone):   p := NULL:   i := i_debut:   while i <> i_fin do     p := p, polygone[i]:     i := voisin_sommet(n, i, 1)   end do:   p := [p, polygone[i_fin]] end proc: #------------------------------------------------------------------- trianguler_polygone_recursif := proc(polygone)   local n, j0, j1, j2, P0, P1, P2, j, polygone_1, polygone_2:   global liste_triangles:   n  := numelems(polygone):   j0 := sommet_gauche(polygone):   j1 := voisin_sommet(n, j0, +1):   j2 := voisin_sommet(n, j0, -1):   P0 := polygone[j0]:   P1 := polygone[j1]:   P2 := polygone[j2]:   j  := sommet_distance_maximale(polygone, P0, P1, P2, [j0, j1, j2]):   if `not`(j::posint) then     liste_triangles := liste_triangles, [P0, P1, P2]:     polygone_1      := nouveau_polygone(polygone,j1,j2):     if numelems(polygone_1) = 3 then       liste_triangles := liste_triangles, polygone_1:     else       thisproc(polygone_1)     end if:   else     polygone_1 := nouveau_polygone(polygone, j0, j ):     polygone_2 := nouveau_polygone(polygone, j , j0):     if numelems(polygone_1) = 3 then       liste_triangles := liste_triangles, polygone_1:     else       thisproc(polygone_1)     end if:     if numelems(polygone_2) = 3 then       liste_triangles := liste_triangles, polygone_2:     else       thisproc(polygone_2)     end if:   end if:   return [liste_triangles]: end proc: #------------------------------------------------------------------- trianguler_polygone := proc(polygone)   trianguler_polygone_recursif(polygone):   return liste_triangles: end proc: #------------------------------------------------------------------- hachurer := proc(shapes, hatch_angle::list, hatch_number::list, hatch_color::list) local A, La, Lp; local N, P, _sides, L_sides, Xshape, ch, rel, p_rel, n, sol, p_range: local AllHatches, window, p, _segment: local NT, ka, N_Hatches, p_range_t, nt, shape, p_hatches, WhatHatches: #----------------------------------------------------------------- # internal functions: # # La(x, y, alpha, p) is the implicit equation of a straight line of angle alpha relatively #                    to the horizontal axis and intercept p # # Lp(x, y, P) is the implicit equation of a straight line passing through points P[1] and P[2] # # interieur_triangle(triangle, M) La := (x, y, alpha, p) -> cos(alpha)*x - sin(alpha)*y + p; Lp := proc(x, y, P::list) (x-P[1][1])*(P[2][2]-P[1][2]) - (y-P[1][2])*(P[2][1]-P[1][1] = 0) end proc; p_range    := [+infinity, -infinity]: NT         := numelems(shapes): AllHatches := NULL: for ka from 1 to numelems(hatch_angle) do   A         := hatch_angle[ka]:   N_Hatches := hatch_number[ka]:   p_range_t := NULL:   _sides    := []:   L_sides   := []:   rel       := []:   for nt from 1 to NT do     shape := shapes[nt]:     # _sides  : two points description of the sides of the shape     # L_sides : implicit equations of the straight lines that support the sides     N        := [1, 2, 3];     P        := [2, 3, 1];     _sides   := [ _sides[] , [ seq([shape[n], shape[P[n]]], n in N) ] ];     L_sides  := [ L_sides[], Lp~(x, y, _sides[-1]) ];     # Inequalities that define the interior of the shape     rel := [ rel[], trianguler:-interieur_forme(shape, [x, y]) ];        # Given the orientation of the hatches we search here the extreme values of     # the intercept p for wich a straight line of equation La(x, y, alpha, p)     # cuts the shape.          p_rel := NULL:          for n from 1 to numelems(L_sides[-1]) do       sol   := solve({La(x, y, A, q), lhs(L_sides[-1][n])} union rel[-1], [x, y]);       p_rel := p_rel, `if`(sol <> [], [rhs(op(1, %)), rhs(op(3, %))], [+infinity, -infinity]);     end do:     p_range_t := p_range_t, evalf(min(op~(1, [p_rel]))..max(op~(2, [p_rel])));     p_range   := evalf(min(op~(1, [p_rel]), op(1, p_range))..max(op~(2, [p_rel]), op(2, p_range)));   end do: # end of the loop over triangles   p_range_t := [p_range_t]:   p_hatches := [seq(p_range, (op(2, p_range)-op(1, p_range))/N_Hatches)]:   # Building of the hatches   #   # This construction is far from being optimal.   # Here again the main goal was to obtain the hatches with a minimum effort   # if algorithmic development.   window      := min(op~(1..shape))..max(op~(1..shape)):   WhatHatches := map(v -> map(u -> if verify(u, v, 'interval'('closed') ) then u end if, p_hatches), p_range_t):   for nt from 1 to NT do     for p in WhatHatches[nt] do       _segment := []:       for n from 1 to numelems(L_sides[nt]) do          _segment := _segment, evalf( solve({La(x, y, A, p), lhs(L_sides[nt][n])} union rel[nt], [x, y]) );       end do;       map(u -> u[], [_segment]);       AllHatches := AllHatches, plot(map(u -> rhs~(u), %), color=hatch_color[ka]):     end do:   end do; end do: # end of loop over hatch angles plots:-display(   PLOT(POLYGONS(polygone, COLOR(RGB, 1\$3))),   AllHatches,   scaling=constrained ) end proc: end module:

Legrand's example (see reference above)

 > global liste_triangles: liste_triangles := NULL:
 > polygone := [[0,0],[0.5,-1],[1.5,-0.2],[2,-0.5],[2,0],[1.5,1],[0.3,0],[0.5,1]]: trianguler:-trianguler_polygone(polygone): PLOT(seq(POLYGONS(u, COLOR(RGB, rand()/10^12, rand()/10^12, rand()/10^12)), u in liste_triangles), VIEW(0..2, -2..2))
 > trianguler:-hachurer([liste_triangles], [-Pi/4, Pi/4], [40, 40], [red, blue])
 > F := (P, a, b) -> map(p -> [p[1]+a, p[2]+b], P):
 > MOTIF  := [[0, 0], [0.05, 0], [0.05, 0.05], [0, 0.05]]; motifs := [ seq(seq(F(MOTIF, 0+i*0.075, 0+j*0.075), i=0..26), j=-14..13) ]: plots:-display(   plot([polygone[], polygone[1]], color=red, filled),   map(u -> plot([u[], u[1]], color=blue, filled, scaling=constrained), motifs) ): texture    := NULL: rel_motifs := map(u -> trianguler:-interieur_forme(u, [x, y]), motifs):    for ref in liste_triangles do   ref;   #   # the three lines below are used to define REF counter clockwise   #   g           := trianguler:-sommet_gauche(ref):   bas         := sort(op~(2, ref), output=permutation);   REF         := ref[[g, op(map(u -> if u<>g then u end if, bas))]];   rel_ref     := trianguler:-interieur_forme(REF, [x, y]): #print(ref, REF, rel_ref);   texture_ref := map(u -> plots:-inequal(rel_ref union u, x=0..2, y=-1..1, color=blue, 'nolines'), rel_motifs):   texture     := texture, texture_ref: end do: plots:-display(   plot([polygone[], polygone[1]], color=red, scaling=constrained),   texture )
 > MOTIF  := [[0, 0], [0.05, 0], [0.05, 0.05], [0, 0.05]]; motifs := [ seq(seq(F(MOTIF, piecewise(j::odd, 0.05, 0.1)+i*0.1, 0+j*0.05), i=-0.2..20), j=-20..20) ]: plots:-display(   plot([polygone[], polygone[1]], color=red, filled),   map(u -> plot([u[], u[1]], color=blue, filled, scaling=constrained), motifs) ): texture    := NULL: rel_motifs := map(u -> trianguler:-interieur_forme(u, [x, y]), motifs):    for ref in liste_triangles do   ref;   g := trianguler:-sommet_gauche(ref):   bas := sort(op~(2, ref), output=permutation);   REF := ref[[g, op(map(u -> if u<>g then u end if, bas))]];   rel_ref     := trianguler:-interieur_forme(REF, [x, y]): #print(ref, REF, rel_ref);   texture_ref := map(u -> plots:-inequal(rel_ref union u, x=0..2, y=-1..1, color=blue, 'nolines'), rel_motifs):   texture     := texture, texture_ref: end do: plots:-display(   plot([polygone[], polygone[1]], color=red, scaling=constrained),   texture )
 > MOTIF  := [[0, 0], [0.4, 0], [0.4, 0.14], [0, 0.14]]: motifs := [ seq(seq(F(MOTIF, piecewise(j::odd, 0.4, 0.2)+i*0.4, 0+j*0.14), i=-1..4), j=-8..7) ]: plots:-display(   plot([polygone[], polygone[1]], color=red, filled),   map(u -> plot([u[], u[1]], color=blue, filled, scaling=constrained), motifs) ): palettes := ColorTools:-PaletteNames(): ColorTools:-GetPalette("HTML"): couleurs := [SandyBrown, PeachPuff, Peru, Linen, Bisque, Burlywood, Tan, AntiqueWhite,      NavajoWhite, BlanchedAlmond, PapayaWhip, Moccasin, Wheat]: nc   := numelems(couleurs): roll := rand(1..nc): motifs_nb      := numelems(motifs): motifs_couleur := [ seq(cat("HTML ", couleurs[roll()]), n=1..motifs_nb) ]: texture    := NULL: rel_motifs := map(u -> trianguler:-interieur_forme(u, [x, y]), motifs):    for ref in liste_triangles do   ref;   g := trianguler:-sommet_gauche(ref):   bas := sort(op~(2, ref), output=permutation);   REF := ref[[g, op(map(u -> if u<>g then u end if, bas))]];   rel_ref     := trianguler:-interieur_forme(REF, [x, y]): #print(ref, REF, rel_ref);   texture_ref := map(n -> plots:-inequal(rel_ref union rel_motifs[n], x=0..2, y=-1..1, color=motifs_couleur[n], 'nolines'), [\$1..motifs_nb]):   texture     := texture, texture_ref: end do: plots:-display(   plot([polygone[], polygone[1]], color=red, scaling=constrained),   texture )
 > MOTIF  := [ seq(0.1*~[cos(Pi/6+Pi/3*i), sin(Pi/6+Pi/3*i)], i=0..5) ]: motifs := [ seq(seq(F(MOTIF, i*0.2*cos(Pi/6)+piecewise(j::odd, 0, 0.08), j*0.3*sin(Pi/6)), i=0..12), j=-6..6) ]: plots:-display(   plot([polygone[], polygone[1]], color=red, filled),   map(u -> plot([u[], u[1]], color=blue, filled, scaling=constrained), motifs) ): motifs_nb      := numelems(motifs): motifs_couleur := [ seq(`if`(n::even, yellow, brown) , n=1..motifs_nb) ]: texture    := NULL: rel_motifs := map(u -> trianguler:-interieur_forme(u, [x, y]), motifs):    for ref in liste_triangles do   ref;   g := trianguler:-sommet_gauche(ref):   bas := sort(op~(2, ref), output=permutation);   REF := ref[[g, op(map(u -> if u<>g then u end if, bas))]];   rel_ref     := trianguler:-interieur_forme(REF, [x, y]): #print(ref, REF, rel_ref);   texture_ref := map(n -> plots:-inequal(rel_ref union rel_motifs[n], x=0..2, y=-1..1, color=motifs_couleur[n], 'nolines'), [\$1..motifs_nb]):   texture     := texture, texture_ref: end do: plots:-display(   plot([polygone[], polygone[1]], color=red, scaling=constrained),   texture )
 >