[perfs] intersections droite/rectangle

Version imprimable

08/04/2008, 18h21
Nemerle

[perfs] intersections droite/rectangle

Tout est dans le titre: je connais différents algos d'intersection d'un point et d'un rectangle quelconque dans le plan (y a eu un sujet pas longtemps par ailleurs).

Quel algo performant d'après vous permet de trouver l'intersection d'un segment et d'un rectangle quelconque dans le plan? ;)
08/04/2008, 18h51
Zavonen
Je crois qu'il suffit de mettre bout à bout tout deux d'algos présentés ici-même il n'y a pas longtemps.
- Test d'appartenance d'un point à un rectangle.
- Intersection de deux segments (proposition récente de Pseudocode).
On teste l'appartenance de chaque extrémité du segment au rectangle.
OUI-OUI --> le segment tout entier.
NON-NON --> Intersections du segment avec chacun des 4 côtés (on peut trouver 0 points ou 2 éventuellement confondus. Donc rien ou un sous-segment, selon.
OUI-NON --> Intersections du segment avec chacun des 4 côtés (on trouve un point, ou éventuellement deux confondus). Sous segment joignant le sommet inclus au point intersection.
Cela dit, on peut peut-être mieux faire.

Pour reprendre l'idée de Zavonen, on peut faire 4 intersections Ligne / segment (chaque segment représente un coté du rectangle).

Chaque intersection est soit vide, 1 point ou 1 segment:

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
 
// Intersection entre le segment [A,B] et la droite O+t.V
static Intersection intersectionLineSegment(Point A, Point B, Point O, Point V) {
	// cas ou AB et V sont colineaires
	double den = (B.x-A.x)*(V.y)-(B.y-A.y)*(V.x);
	if (den==0)  { 
		// O est sur la droite (A,B) -> intersection = tout le segment
		double det = (B.x-A.x)*(O.y-A.y)-(B.y-A.y)*(O.x-A.x);
		if (det==0) return new IntersectionSegment(A,B);	
 
		// sinon pas d'intersection
		return null;     
	}
 
	// calcul l'intersection entre les 2 droites (A,B) et O+t.V
	double num = (A.y-O.y)*(V.x)-(A.x-O.x)*(V.y);
	if (den<0) {num=-num; den=-den;}
 
	// l'intersection est avant A ou après B -> pas d'intersection
	if (num<0)   return null;
	if (num>den) return null;
 
	// l'intersection est le point A
	if (num==0) return new IntersectionPoint(A); 
 
	// l'intersection est le point B
	if (num==den) return new IntersectionPoint(B); 
 
	/* note: on peut ignorer le cas ou l'intersection est le
	 * point B car, dans le cas du rectangle, ce point sera
	 * l'origine d'un autre segment.		
	 */
 
	// sinon l'intersection est entre A et B
	double r = num/den;
	Point p = new Point(A.x+r*(B.x-A.x), A.y+r*(B.y-A.y) );
	return new IntersectionPoint(p);
}

Le résultat final est l'union des 4 intersections.

(je vous laisse écrire le code :aie:)

08/04/2008, 20h18
souviron34

[un petit peu HS]

ca me fait penser, je ne sais pas si vous connaissez cette excellente et simplissime routine pour effectuer des calculs de contours , CONREC, de P.Burke (http://local.wasp.uwa.edu.au/~pbourke/papers/conrec/)...

Explication simple, application egalement...

Si il y a des valeurs attribuees aux sommets du rectangle et aux extremtites du segment, application pour le cas present.

Mais sinon, j'ai ete tellement ebloui par la beaute et simplicite de cette methode que je voulais la faire partager...

[/un petit peu HS]
09/04/2008, 00h04
pseudocode

[un petit peu HS]une variante du marching square ?[/un petit peu HS]
09/04/2008, 08h30
Zavonen

On peut aussi raisonner ainsi.
Commencer par traiter le problème de l'intersection d'une droite avec un rectangle. Cela a déjà été fait ici il n'y a pas longtemps. Le critère d'intersection est simple si on compare la distance du centre de gravité du rectangle avec la demi-diagonale.
On peut donc déjà faire une première discrimination presque immédiate:
La droite (AB) coupe-t-elle le rectangle ?
Si NON problème terminé.
Si OUI determiner le segment intersection (problème à priori plus simple que le problème initial, cela peut se faire juste avec ce bon vieux Pythagore).
Ensuite il ne reste plus qu'à déterminer l'intersection de deux segments d'une même droite, ce qui revient à des comparaisons de coordonnées.
09/04/2008, 10h01
Nemerle

Merci pour ces retours.

Je précise (juste pour info) que le rectangle est plein: c'est donc un sous-segment du segment initial que l'on recherche.

Autre idée possible: pourquoi ne pas changer de repère pour "remettre" le rectangle parallèle aux axes? L'intersection est alors triviale à définir...
09/04/2008, 10h14
Zavonen

Citation:

Autre idée possible: pourquoi ne pas changer de repère pour "remettre" le rectangle parallèle aux axes? L'intersection est alors triviale à définir...

C'était une évidence (sous-entendue). Maintenant tout dépend de ce qu'on entend par 'trivial', ce qui est 'trivial' pour l'un ne l'est pas pour l'autre (enfin ça aussi c'est une trivialité).
09/04/2008, 11h06
pseudocode

Citation:

Envoyé par Nemerle

Je précise (juste pour info) que le rectangle est plein: c'est donc un sous-segment du segment initial que l'on recherche.

Ah, si tu nous ne dis pas tout. Dans ce cas c'est un poil plus facile, surtout avec l'equation parametrique de la droite.

Il suffit de conserver le "tmin" et "tmax" des 4 intersections ligne/segment: l'intersection finale est l'ensemble des points O+t.V avec t dans [tmin,tmax].

Citation:

Autre idée possible: pourquoi ne pas changer de repère pour "remettre" le rectangle parallèle aux axes? L'intersection est alors triviale à définir...

En théorie: oui. En pratique: les erreurs de calculs sur les flottants faussent les résultats.
09/04/2008, 11h29
Zavonen
Citation:

En théorie: oui. En pratique: les erreurs de calculs sur les flottants faussent les résultats.

OUI, à tous les niveaux, même pour une banale intersection de deux droites. Il faut donc être convaincu une fois pour toutes que les brillantes idées de Descartes jointes au génie de Von Neumann n'offre que des approximations de problèmes extrêmement simples ou autrefois la règle et le compas faisaient merveille.
Pour mettre tout au clair:
- Si le repère est centré sur le centre de gravité du rectangle.
- Si les vecteurs unitaires sont parallèles aux côtés du rectangle.
- Si L est la longueur et l la largeur du rectangle.
La distance d de O à la droite d'équation ax+by+c=0 est
abs(c)/sqrt(a*a+b*b)
La demi diagonale e est sqrt(L*L+l*l)/2
Donc si d>e pas de solution.
Si d<e on détermine les points d'intersection en examinant d'abord les cas a=0 et b=0 (droites parallèles aux côtés du triangle) les solutions sont alors évidentes.
Sinon il suffit de remplacer dans l'équation de la droite x et y par -L/2 +L/2, -l/2, +l/2 et de calculer l'autre coordonnée pour avoir les points d'intersection (en supposant que c'est la longueur qui détermine l'axe des abscisses).
Donc en résumé:
Si on donne deux points A et B, on commence par déterminer l'équation de la droite (AB) en annulant par exemple det(MA,MB).
On détermine par la méthode ci dessus ,s'ils existent, les deux points d'intersection (éventuellement confondus) I et J du périmètre du rectangle avec (AB).
La réponse est l'intersection des segments [AB] et [IJ] que l'on détermine par une simple comparaison des coordonnées des points.
Faire plus simple, je ne vois pas comment.

Citation:

Envoyé par pseudocode

double den = (B.x-A.x)*(V.y)-(B.y-A.y)*(V.x);
if (den==0) {
// O est sur la droite (A,B) -> intersection = tout le segment

Au fait, CA, c'est fô...

Maintenant, ma soluce (crade :mouarf:, sous excel) proche de celle de Pseudocode (segment (x1,y1)--->(x2,y2), et les cotés du rectangle sont donnés par o(i,j)).

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
nbsol = 0
 
For i = 0 To 3
    ok = False
    u(1) = x1 - o(i, 1)
    u(2) = y1 - o(i, 2)
    v(1) = x2 - o(i, 1)
    v(2) = y2 - o(i, 2)
    w(1) = o(i + 1, 1) - o(i, 1)
    w(2) = o(i + 1, 2) - o(i, 2)
 
    den = u(1) * v(2) - v(2) * u(1)
    det = w(1) * u(2) - w(2) * u(1)
 
   If den = 0 And det = 0 And det <> den Then '====== (O(i)O(i+1)) = (AB)
        ta = u(1) / w(1)
        tb = v(1) / w(1)
        If Not ((ta < 0 Or ta >= 1) And (tb < 0 Or tb >= 1)) Then
            If ta > tb Then
                tc = tb
                tb = ta
                ta = tc
            End If
            If ta < 0 Then ta = 0
            If tb > 1 Then tb = 1
            sol(1) = ta
            sol(2) = tb
            nbsol = 2
            i = 4
        End If
    Else
        det1 = w(1) * ab(2) - w(2) * ab(1)
        If det1 <> 0 Then
            mu = -det / det1
        Else
            mu = -1
        End If
        sol(nbsol) = mu
        nbsol = nbsol + 1
    End If
Next i
For i = 0 To nbsol - 1
    If (sol(i) >= 0 And sol(i) <= 1) Then
        solgood(nbgood) = sol(i)
        nbgood = nbgood + 1
    End If
Next i
If nbgood = 1 Then
    If isinside(x1, y1, lx, ly, rx, ry, tx, ty, bx, by) Then
        solgood(nbgood) = 0
        nbgood = nbgood + 1
    End If
    If isinside(x2, y2, lx, ly, rx, ry, tx, ty, bx, by) Then
        solgood(nbgood) = 1
        nbgood = nbgood + 1
    End If
End If
If nbgood = 2 Then
    intersect = Abs(solgood(1) - solgood(0))
Else
    intersect = 0
End If

Citation:

Envoyé par Nemerle

Au fait, CA, c'est fô...

:koi: ?
Code:

1 2 3 // cas ou AB et V sont colineaires // O est sur la droite (A,B) -> intersection = tout le segment

09/04/2008, 14h55
Nemerle
Citation:
Envoyé par pseudocode

:koi: ?

Code:

1 2 3 // cas ou AB et V sont colineaires // O est sur la droite (A,B) -> intersection = tout le segment
:koi: ?:koi: ? ;)

Tu calcules le determinant des vecteurs AB et V: ça te dis juste effectivement que AB et V sont colinéaires, mais pas que 0 est sur (AB) !

Nàn?
09/04/2008, 15h08
pseudocode

Citation:

Envoyé par Nemerle

Tu calcules le determinant des vecteurs AB et V: ça te dis juste effectivement que AB et V sont colinéaires, mais pas que 0 est sur (AB) !

Nàn?

Arf. C'est juste que mes lignes de commentaires sont AVANT les calculs, et non pas après. :mrgreen:

Voici la traduction en code Python des idées que j'ai développées.
On raisonne ici avec des points exprimés en pixels, il est donc possible de lister de façon exhaustive les points d'un segment.
Il reste peut être (certainement) des bugs mais :
On ne change pas de repère
On n'utilise que 4 opérations

Code:

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
#objet point
class Point:
    #initialisation à partir des coordonnées
    def __init__(self, a,b):
        self.x=a
        self.y=b
    #translaté d'un point par un vecteur
    def __add__(self,V):
        return Point(self.x+V.x,self.y+V.y)
    #vecteur comme différence de deux points
    def __sub__ (self,other):
        return Vecteur(other,self)
 
#objet vecteur
class Vecteur:
    #initialisation à partir de deux points
    def __init__(self,A,B):
        self.x=B.x-A.x
        self.y=B.y-A.y
    #somme de deux vecteurs
    def __add__(self,other):
        X=Point(self.x,self.y)
        return X+other
    #produit saclaire
    def __mul__(self,other):
        return self.x*other.x+self.y*other.y
    #norme
    def __abs__(self):
        return self.x*self.x+self.y*self.y
    #produit d'un vecteur par un scalaire
    def scal(self,m):
        x=int(self.x*m)
        y=int(self.y*m)
        return Point(x,y)-Point(0,0)
 
#objet rectangle
class Rectangle:
    #constructeur avec 3 sommets
    def __init__(self,X,Y,Z):
        self.A=X
        self.B=Y
        self.D=Z
        self.C= (X+(Y-X))+(Z-X)
 
#objet segment
class Segment:
    #initialisation avec les point extrêmes
    def __init__(self,P,Q):
        self.M=P
        self.N=Q
    # le point de paramètre m m dans [0,1]
    def Param(self,m):
        V=(self.N-self.M).scal(m)
        return self.M+V
 
# coordonnées de V dans la base (V1,V2) les deux vecteurs étant supposés orthogonaux
def coordonnees (V,V1,V2):
    return (V*V1)/(V1*V1),(V*V2)/(V2*V2)
 
#teste si P est dans le rectangle R
def IsIn(P,R):
    V1=R.B-R.A
    V2=R.D-R.A
    V= P-R.A
    return 0<=coordonnees(V,V1,V2)[0]<=1 and 0<=coordonnees(V,V1,V2)[1]<=1
 
#imprime tous les points dans R inter S
def PrintInter(R,S):
    V=S.N-S.M
    if V.x > V.y : Max=V.x
    else : Max=V.y
    for i in range(0,Max):
        m= i/float(Max)
        Mm=S.Param(m)
        if IsIn(Mm,R):
            print Mm.x,Mm.y
 
def main():
    A=Point(0,100)
    B=Point(100,0)
    D=Point(200,300)
    R=Rectangle(A,B,D)
    C=R.C
    T1=Point(150,300)
    T2=Point(200,0)
    S=Segment(T1,T2)
    PrintInter(R,S)
 
 
 
if __name__ == '__main__':
    main()

13/04/2008, 11h30
Alikendarfen

Je n'ai pas compris pourquoi le changement de repère fausse les calculs ?
Pourquoi le ferait il plus que tout autre calcul ?

Citation:

Autre idée possible: pourquoi ne pas changer de repère pour "remettre" le rectangle parallèle aux axes? L'intersection est alors triviale à définir...

Il faudrait plutôt mettre le segment dans le repère (à choisir) du rectangle ?
13/04/2008, 12h23
Zavonen

Citation:

Envoyé par Alikendarfen

Je n'ai pas compris pourquoi le changement de repère fausse les calculs ?
Pourquoi le ferait il plus que tout autre calcul ?
?

Oui c'est ce que j'ai déjà fait remarquer TOUS les calculs sont faux:

Citation:

OUI, à tous les niveaux, même pour une banale intersection de deux droites. Il faut donc être convaincu une fois pour toutes que les brillantes idées de Descartes jointes au génie de Von Neumann n'offre que des approximations de problèmes extrêmement simples ou autrefois la règle et le compas faisaient merveille.