Projection de vecteurs sur axes.

Invité · 29/06/2014, 21h26

Salut, je cherche à déterminer si un point est dans un triangles, pour cela je suis entrain de mettre au poins un algorithme généraliste de tels manière qu'il me serve aussi pour déterminer les intersections entre polygones et rayon/polygone.

Pour cela j'essaye de tout ramener à une collision de type cercle/rayon, cercle/cercle ou bien cercle/point à l'aide de projections.

Cependant il y a quelque chose que je ne comprend pas, je n'ai pas les bons résultats pour les projections, voici une image qui montre ce que je cherche à faire :
Nom : img1.jpg
Affichages : 159
Taille : 35,0 Ko

Nom : img1.jpg
Affichages : 159
Taille : 35,0 Ko

Pour la première partie qui se charge d'appliquer l'algorithme de voronoi pour trouver dans quel région se trouve le point p, ça pas de problème :

Code c++ :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
 
Vec3f v = axis - center;
bool found = false;
Vec3f v1, v2, n;
for (int i = 0; i < vertices.size() && !found; i++) {
      v1 = vertices[i] - center;
      v2 = vertices[(i+1 == vertices.size()) ? 0 : i + 1] - center;
      float angle1 = v.getAngleBetween(v1);
      float angle2 = v.getAngleBetween(v2);
      n = normals[i] - center;
     float dp = n.dot(v);
     if (dp >= 0 && (angle1 <= 0 && angle2 > 0 || angle1 > 0 && angle2 <= 0))
         found = true;
    }
}

En fait le problème survient plus loin, lorsque j'essaye de projeter le vecteur v2 sur la normale (je dirais la bissectrice plutôt dans mon cas) n1.

Supposons que j'ai un triangle de sommets (0, 0, 0) ; (50, 100, 0) et (-50, 100, 0).

Le problème est que lorsque je veux projeter v2 sur la bissectrice n1 comme ceci : (v2 et n1 sont normalisé pour avoir le cosinus de l'angle entre n1 et v2)

Code c++ :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
 
float projection = n1.produitScalaire(v2) * v2.longueur();

Je me retrouve avec une projection de 36.98 avec mon exemple.

Et cette valeur est différente de la longueur de n1 qui est de 30.0463.

Hors que ces 2 valeurs devraient être exactement les mêmes car le triangle formé par le centre, par v2 et par l'intersection de la bissectrice avec le côté du triangle est un triangle rectangle.

Bref je ne comprends pas du tout pourquoi ces 2 valeurs ne sont pas les mêmes.

Sinon par la suite je recherche l'intersection à l'aide du théorème de thalès :

Code c++ :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
 
float max = v.magnitude();
float p = v.projOnAxis(n);
float ratio = n.magnitude() / p;
float min = max - v.magnitude() * ratio;

Du coup quand je veux tester cela avec le point (0, 0, 0) il m'indique que ce point n'est pas dans le triangle, hors qu'il l'est vu que c'est un des 3 sommets.

La classe polygone :

Code c++ :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
 
bool BoundingPolyhedron::isPointInside(Vec3f p) {
    float c = (p - center).magnitude();
    Vec2f pr = Computer::projectShapeOnAxisFromCenter(p, center, points, normals);
    float r = pr.y - pr.x;
    std::cout<<"r : "<<r<<" c : "<<c<<std::endl;
    if (r - c < 0)
        return false;    //}
    return true;
}

Je vois pas ou est mon erreur. :/

Pourquoi est ce que la projection de v2 sur n1 ne me donne pas la longueur de n1 comme cela le devrait ???

Sinon, voici le code source complet :

Code c++ :

Sélectionner tout - Visualiser dans une fenêtre à part

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
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
 
#include "../../../include/odfaeg/Physics/boundingPolyhedron.h"
#include "../../../include/odfaeg/Physics/boundingSphere.h"
#include "../../../include/odfaeg/Physics/boundingEllipsoid.h"
#include "../../../include/odfaeg/Physics/boundingBox.h"
#include "../../../include/odfaeg/Physics/orientedBoundingBox.h"
namespace odfaeg {
BoundingPolyhedron::BoundingPolyhedron () {
}
void BoundingPolyhedron::addPoint(Vec3f point) {
    bool contains = false;
    for (unsigned int i = 0; i < points.size() && !contains; i++) {
        if (points[i] == point)
            contains = true;
    }
    if (!contains) {
        points.push_back(point);
        computeCenter();
        computeNormals();
    }
}
bool BoundingPolyhedron::intersects (BoundingSphere &bs) {
    return bs.intersects (*this);
}
bool BoundingPolyhedron::intersects (BoundingEllipsoid &be) {
    return be.intersects (*this);
}
bool BoundingPolyhedron::intersects (BoundingBox& bx) {
    return bx.intersects (*this);
}
bool BoundingPolyhedron::intersects (OrientedBoundingBox &obx) {
    return obx.intersects (*this);
}
bool BoundingPolyhedron::intersects (BoundingPolyhedron &bp) {
     std::vector<Vec3f> axis1 = normals;
     std::vector<Vec3f> axis2 = bp.normals;
     for (unsigned int i = 0; i < normals.size(); i++) {
        Vec2f p1 = Computer::projectShapeOnAxis(normals[i], points);
        Vec2f p2 = Computer::projectShapeOnAxis(normals[i], bp.points);
        if (!Computer::overlap(p1, p2))
            return false;
     }
     for (unsigned int i = 0; i < bp.normals.size(); i++) {
        Vec2f p1 = Computer::projectShapeOnAxis(bp.normals[i], points);
        Vec2f p2 = Computer::projectShapeOnAxis(bp.normals[i], bp.points);
        if (!Computer::overlap(p1, p2))
            return false;
     }
     return true;
}
bool BoundingPolyhedron::isPointInside(Vec3f p) {
    float c = (p - center).magnitude();
    Vec2f pr = Computer::projectShapeOnAxisFromCenter(p, center, points, normals);
    float r = pr.y - pr.x;
    std::cout<<"r : "<<r<<" c : "<<c<<std::endl;
    if (r - c < 0)
        return false;    //}
    return true;
}
int BoundingPolyhedron::nbPoints () {
    return points.size();
}
void BoundingPolyhedron::computeCenter () {
    Vec3f sum(0.f, 0.f, 0.f);
    for (unsigned int i = 0; i < points.size(); i++) {
        sum += points[i];
    }
    center = sum / points.size();
    std::array<std::array<float, 2>, 3> store = Computer::getExtends(points);
    size = Vec3f(store[0][1] - store[0][0], store[1][1] - store[1][0], store[2][1] - store[2][0]);
}
void BoundingPolyhedron::computeNormals() {
    normals.clear();
    for (unsigned int i = 0; i < points.size(); i++) {
        Vec3f p1 = points[i];
        Vec3f p2 = (i + 1 == points.size()) ? points[0] : points[i+1];
        Vec3f n = (p1 + p2) * 0.5f;
 
        /*Vec3f n (p1.y - p2.y, p2.x - p1.x, 0);
        if (n.dot(center- p1) > 0)
            n = -n;*/
        normals.push_back(n);
 
    }
}
std::vector<Vec3f> BoundingPolyhedron::getNormals() {
    return normals;
}
std::vector<Vec3f> BoundingPolyhedron::getPoints() {
    return points;
}
BoundingPolyhedron::~BoundingPolyhedron () {
    normals.clear();
    points.clear();
}
bool BoundingPolyhedron::operator== (const BoundingPolyhedron &bp) {
    if(points.size() != bp.points.size())
        return false;
    for (unsigned int i = 0; i < points.size(); i++) {
        if (!(points[i] == bp.points[i])) {
            return false;
        }
    }
    return true;
}
Vec3f BoundingPolyhedron::getPoint(unsigned  int index) {
    if (index >= 0 && index < points.size())
        return points[index];
    else
        return Vec3f (0, 0, 0);
}
 
Vec3f BoundingPolyhedron::getCenter () {
    return center;
}
Vec3f BoundingPolyhedron::getSize() {
    return size;
}
Vec3f BoundingPolyhedron::getPosition() {
    return Vec3f (center.x - size.x * 0.5f, center.y - size.y * 0.5f, center.z - size.z * 0.5f);
}
}

La classe vec4f :

Code c++ :

Sélectionner tout - Visualiser dans une fenêtre à part

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
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
 
#include "../../../include/odfaeg/Math/vec4.h"
#include "../../../include/odfaeg/Math/matrix4.h"
namespace odfaeg {
using namespace std;
using namespace sf;
const Vec3f Vec3f::xAxis(1.f, 0.f, 0.f);
const Vec3f Vec3f::yAxis(0.f, 1.f, 0.f);
const Vec3f Vec3f::zAxis(0.f, 0.f, 1.f);
 
ostream& operator<< (ostream &out, const Vec3f &vec3) {
    out<<"x : "<<vec3.x<<" y : "<<vec3.y<<" z : "<<vec3.z<<" w : "<<vec3.w<<endl;
    return out;
}
Vec3f::Vec3f () : Vector3f (0.f, 0.f, 0.f) {
    w = 1.f;
}
 
Vec3f::Vec3f(float x, float y, float z) : Vector3f (x, y, z) {
    w = 1.f;
}
Vec3f::Vec3f (float x, float y, float z, float w) : Vector3f(x, y, z) {
    this->w = w;
}
 
void Vec3f::set (float x, float y, float z) {
    this->x = x;
    this->y = y;
    this->z = z;
}
 
Vec3f Vec3f::normalizeToVec3 () const {
 
    Vector3f vec3 = getVec3sf();
    if (w != 0.f) {
        return Vec3f (vec3.x / w, vec3.y / w, vec3.z / w);
 
    }
    else
        return *this;
}
bool Vec3f::isNulVector () const {
    return x==0 && y==0 && z==0;
}
Vec3f Vec3f::operator- () const {
    return Vec3f (-this->x, -this->y, -this->z);
}
Vec3f Vec3f::operator+ (const Vec3f &other) {
 
    Vector3f vec3 = this->getVec3sf() + other.getVec3sf();
    return Vec3f (x + other.x, y + other.y, z + other.z);
 
}
Vec3f Vec3f::operator- (const Vec3f &other) {
    sf::Vector3f v1 = this->getVec3sf();
    sf::Vector3f v2 = other.getVec3sf();
    sf::Vector3f vec3 = v1 - v2;
    return Vec3f (x - other.x, y - other.y, z - other.z);
}
Vec3f Vec3f::operator* (const Vec3f &other) {
    float x = this->x * other.x;
    float y = this->y * other.y;
    float z = this->z * other.z;
    return Vec3f (x, y, z);
}
Vec3f Vec3f::operator* (const float scale) {
    sf::Vector3f v1 = this->getVec3sf();
    sf::Vector3f vec3 = v1 * scale;
    return Vec3f (x * scale, y * scale, z * scale);
}
Vec3f Vec3f::operator/ (const float scale) {
    sf::Vector3f v1 = this->getVec3sf();
    sf::Vector3f vec3 = v1 / scale;
    return Vec3f (vec3.x, vec3.y, vec3.z);
}
Vec3f Vec3f::operator/ (const Vec3f &other) {
    try {
        if (other.x == 0 || other.y == 0 || other.z == 0)
            throw Erreur(0, "Division by zero.", 1);
        else
            return Vec3f (this->x / other.x, this->y / other.y, this->z / other.z);
    } catch (Erreur err) {
        cerr<<err.what()<<endl;
    }
}
Vec3f& Vec3f::operator+= (const Vec3f &other) {
 
    this->x += other.x;
    this->y += other.y;
    this->z += other.z;
    return *this;
 
}
void Vec3f::operator-= (const Vec3f &other) {
   this->x -= other.x;
   this->y -= other.y;
   this->z -= other.z;
}
void Vec3f::operator*= (const Vec3f &other) {
   this->x *= other.x;
   this->y *= other.y;
   this->z *= other.z;
}
void Vec3f::operator*= (const float scale) {
   this->x *= scale;
   this->y *= scale;
   this->z *= scale;
}
void Vec3f::operator/= (const Vec3f &other) {
    try {
        if (other.x == 0 || other.y == 0 || other.z == 0)
            throw Erreur(0, "Division by zero.", 1);
        else {
            this->x /= other.x;
            this->y /= other.y;
            this->z /= other.z;
        }
    } catch (Erreur err) {
        cerr<<err.what()<<endl;
    }
}
 
 
Vec3f& Vec3f::operator= (const Vec3f &other) {
    x = other.x;
    y = other.y;
    z = other.z;
    w = other.w;
    return *this;
}
bool Vec3f::operator== (const Vec3f &other) {
    return x==other.x && y==other.y && z==other.z;
}
 
float Vec3f::computeDist (const Vec3f &other) {
    return Math::sqrt(Math::power(x - other.x, 2) + Math::power(y - other.y, 2) + Math::power(z - other.z, 2));
}
float Vec3f::magnitude () const {
    return Math::sqrt(Math::power(x, 2) + Math::power(y, 2) + Math::power(z, 2));
}
float Vec3f::magnSquared () const {
    return Math::power(x, 2) + Math::power(y, 2) + Math::power(z, 2);
}
Vec3f Vec3f::normalize () const {
    float length = magnitude();
 
    Vec3f v3f;
    if (length==0) {
        v3f.set(0.f, 0.f, 0.f);
    } else {
        float x = this->x / length;
        float y = this->y / length;
        float z = this->z / length;
        v3f.set (x, y, z);
    }
 
    return v3f;
}
float Vec3f::dot2 (const Vec3f &other) {
    if (isNulVector ())
        return 0.f;
    return x * other.x + y * other.y + z * other.z;
 
}
float Vec3f::dot3(const Vec3f &other) {
    Vec3f v (x + other.x, y + other.y, z + other.z);
    return 0.5f * Math::power(v.magnitude(), 2) - Math::power(magnitude(), 2) - Math::power(other.magnitude(), 2);
}
float Vec3f::dot (const Vec3f &other) const {
    if(isNulVector ())
        return 0;
   Vec3f n1 = normalize();
   Vec3f n2 = other.normalize();
   return (n1.x * n2.x + n1.y * n2.y + n1.z * n2.z);
}
float Vec3f::getAngleBetween (const Vec3f &other) {
    if(isNulVector() || other.isNulVector())
        return 0;
    float dotProduct = dot(other);
    Vec3f n = cross(other);
    Matrix4f matrix (x, y, z, 0, other.x, other.y, other.z, 0, n.x, n.y, n.z, 0, 0, 0, 0, 1);
    float fDet = matrix.getDet();
    if (fDet >= 0)
        return Math::acosinus(dotProduct);
    else
        return -Math::acosinus(dotProduct);
}
Vec3f Vec3f::cross (const Vec3f &other) {
 
    return Vec3f (y * other.z - z * other.y, z * other.x - x * other.z, x * other.y - y * other.x);
}
float Vec3f::projOnAxis (const Vec3f &other) {
    float dp = dot(other);
    return dp * magnitude();
}
Vec3f Vec3f::projOnVector(Vec3f other) {
    Vec3f proj;
    float dp = dot2(other);
    proj.x = dp / other.magnSquared() * other.x;
    proj.y = dp / other.magnSquared() * other.y;
    proj.z = dp / other.magnSquared() * other.z;
    return proj;
}
bool Vec3f::isOpposite (const Vec3f &other) const {
    if (isNulVector() || other.isNulVector())
        return false;
    if ((x > 0.f && other.x < 0.f) || (x < 0.f && other.x > 0.f)
        && (y > 0.f && other.y < 0.f) || (y < 0.f && other.y > 0.f)
        && (z > 0.f && other.z < 0.f) || (z < 0.f && other.z > 0.f))
        return true;
    return false;
}
float* Vec3f::getVec3 () const {
    float v[3] = {x, y, z};
    return v;
}
Vector3f Vec3f::getVec3sf () const {
    return Vector3f(x, y, z);
}
 
Vec3f Vec3f::fabs () const {
    float x = Math::abs(this->x);
    float y = Math::abs(this->y);
    float z = Math::abs(this->z);
    return Vec3f (x, y, z);
}
float& Vec3f::operator[](int i) {
    if (i == 0)
        return x;
    else if (i == 1)
        return y;
    else if (i == 2)
        return z;
    else
        return w;
}
}

La classe math :

Code c++ :

Sélectionner tout - Visualiser dans une fenêtre à part

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
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
 
#include "../../../include/odfaeg/Math/maths.h"
#include "../../../include/odfaeg/Math/vec4.h"
/*Fonctions élémentaires.*/
//Donne l'arc cosinus d'une valeur.
 
namespace odfaeg {
std::mt19937 Math::mrs = std::mt19937(static_cast<unsigned long>(std::time(nullptr)));
float Math::random(float min, float max) {
 
    std::uniform_real_distribution<float> distribution(min, max);
    return distribution(mrs);
}
 
float Math::acosinus(float value) {
   float result;
   if(-1.f < value) {
    if(value < 1.f)
       result = (float) acos(value);
    else
       result = 0.f;
   } else
    result = PI;
   return (result < EPSILON && result > -EPSILON) ? 0 : result;
}
//Donne l'arc sinus d'une valeur.
float Math::asinus(float value) {
    float result;
    if(-1.f < value) {
        if(value < 1.f)
            result = (float) asin(value);
        else
            result = HALF_PI;
    } else
        result = -HALF_PI;
    return (result < EPSILON && result > -EPSILON) ? 0 : result;
}
//Donne l'arc tangeante d'une valeur.
float Math::atang(float value) {
    float result = (float) atan(value);
    return (result < EPSILON && result > -EPSILON) ? 0 : result;
}
//Donne le sinus d'un angle donné en radian.
float Math::sinus (float value) {
    float result, radians = value / TWO_PI;
    if(abs(radians) > PI)
        radians -= TWO_PI;
    if(abs(radians) > HALF_PI)
        radians = PI - radians;
    if(abs(radians) <= PI / 4) {
        result = (float) sin(value);
    } else {
        result = (float) cos (PI/2-value);
    }
    return (result < EPSILON && result > -EPSILON) ? 0 : result;
}
//Donne le cosinus d'un angle donnée en radian.
float Math::cosinus (float value) {
    float result = (float) sin(value+HALF_PI);
    return (result < EPSILON && result > -EPSILON) ? 0 : result;
}
//Donne la tangeante d'un angle donné en radians.
float Math::tang (float value) {
    float result = (float) tan(value);
    return (result < EPSILON && result > -EPSILON) ? 0 : result;
}
//Donne la valeur absolue d'un nombre.
float Math::abs(float value) {
    if (value >= 0)
        return value;
    return -value;
}
//Renvoie la racine carrée d'un nombre.
double Math::sqrt(float value){
    return std::sqrt(value);
}
//Donne l'inverse de la racine carrée d'un nombre.
float Math::inversSqrt(float value) {
    return 1.f / sqrt(value);
}
//Donne le logarithme d'un nombre. (En base 10.)
float Math::log10 (float value) {
    return (float) log(value);
}
//Donne le logarithme d'un nombre en base base.
float Math::logn (float value, int base) {
    return (float) (log(value)/log(base));
}
//Donne le nombre à la puissance n.
double Math::power (float value, float exp) {
    return (float) pow(value, exp);
}
//Converti un angle en radian.
float Math::toRadians(float value) {
    return value * DEG_TO_RAD;
}
//Convertis un angle en degrer.
float Math::toDegrees(float value) {
    return value * RAD_TO_DEG;
}
//Arrondis un nombre à la précision p.
float Math::round(float value, int p) {
    int mult = (int) pow(10, p + 1);
	int numberToRound = (int) (value * mult);
 
	int lastChiffer = numberToRound % 10;
 
	if (numberToRound > 0) {
		if (lastChiffer >= 5)
			numberToRound += 10;
	} else {
		if (lastChiffer <= -5)
			numberToRound -= 10;
	}
 
	numberToRound = numberToRound - lastChiffer;
	return numberToRound / mult;;
}
//Renvoie l'exponetielle d'un nombre.
float Math::exp (float value) {
    return (float) std::exp(value);
}
//Convertis des coordonnée polaire en coordonée cartésinnes.
Vec3f Math::toCartesian (float teta, float phi) {
    float rTemp = Math::cosinus(phi);
    float y = Math::sinus(phi);
    float z = -rTemp * Math::cosinus(teta);
    float x = -rTemp * Math::sinus(teta);
    return Vec3f (x, y, z);
 
}
 
 
int Math::clamp (int value, int min, int max) {
    if (value < min)
        value = min;
    if (value >= max);
        value = max - 1;
 
    return value;
}
}

Le main :

Code c++ :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
 
int main() {
odfaeg::BoundingPolyhedron t2;
t2.addPoint(odfaeg::Vec3f(0, 0, 0));
t2.addPoint(odfaeg::Vec3f(50, 100, 0));
t2.addPoint(odfaeg::Vec3f(-50, 100, 0));
odfaeg::Vec3f pts(0, 0, 0);
std::cout<<t2.isPointInside(pts);
return 0;
}

Et la classe qui me pose problème :

Code c++ :

Sélectionner tout - Visualiser dans une fenêtre à part

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
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
141
142
143
144
145
146
147
148
149
150
151
152
153
154
155
156
157
158
159
160
161
162
163
164
165
166
167
168
169
170
171
172
173
174
175
176
177
178
179
180
181
182
183
184
185
186
187
188
189
190
191
192
193
194
195
196
197
198
199
200
201
202
203
204
205
206
207
208
209
210
211
212
213
214
215
216
217
218
219
220
221
222
223
224
225
226
227
228
229
230
231
232
233
234
235
236
237
238
239
240
241
242
243
244
245
246
247
248
249
250
251
252
253
254
255
256
257
258
259
260
261
262
263
264
265
266
267
268
269
270
271
272
273
274
275
276
277
278
279
280
281
282
283
284
285
286
287
288
289
290
291
292
293
294
295
296
297
298
299
300
301
302
303
304
305
306
307
308
309
310
311
312
313
314
315
316
317
318
319
320
321
322
323
324
325
326
327
328
 
#ifndef COMPUTER
#define COMPUTER
#include "matrix3.h"
#include "vec4.h"
#include <limits.h>
#include <array>
#include "../Physics/ray.h"
#include <vector>
#include "export.hpp"
/**
  *\namespace odfaeg
  * the namespace of the Opensource Development Framework Adapted for Every Games.
  */
namespace odfaeg {
/**
  * \file computer.h
  * \class Computer
  * \brief Compute the min and max points of a vector or an array of points.
  * \author Duroisin.L
  * \version 1.0
  * \date 1/02/2014
  *
  * Manage a computer to check the min and max x,y coordinates from a set of points or the average.
  * Return the average or a 2 dimentionnal array containing the min and max coordinates.
  * The index of min and max are : 0,0 for the min x; 0, 1 for the max x; 1,0 for the min Y, etc...
  */
class ODFAEG_MATH_API Computer {
    public :
    /**\fn Vec2f getMoy (std::vector<Vec2f*> verts)
    *  \brief get the average from a set of 2D vectors.
    *  \param a list of vectors.
    *  \return the average vector.
    */
    static Vec2f getMoy (std::vector<Vec2f> verts);
    static Vec3f getMoy (std::vector<Vec3f> verts);
    /**\fn std::array<std::array<float, 3>,2> getExtends (std::vector<Vec3f*> verts);
    *  \brief get the minimum and the maximum x, y and z from a set of 3D vectors.
    *  \param the list of the vectors.
    *  \return an array containing the minimum and maximum values.
    */
    static std::array<std::array<float, 2>,3> getExtends (std::vector<Vec3f> verts);
    /**\fn std::array<std::array<float, 2>,2> getExtends (std::vector<Vec2f*> verts);
    *  \brief get the minimum and the maximum x and y from a set of 2D vectors.
    *  \param the list of the vectors.
    *  \return an array containing the minimum and maximum values.
    */
    static std::array<std::array<float, 2>,2> getExtends (std::vector<Vec2f> verts);
    /**\fn std::array<std::array<float, 2>,2> getExtends (const std::arrayr<Vec2f, N>& verts);
    *  \brief get the minimum and the maximum x and y from an array of 2D vectors.
    *  \param the array of the vectors.
    *  \return an array containing the minimum and maximum values.
    */
    template <std::size_t N>
    static std::array<std::array<float, 2>,2> getExtends (const std::array<Vec2f, N>& verts) {
        float minX = 0;
        float maxX = 0;
        float minY = 0;
        float maxY = 0;
        if (verts.size() > 0) {
            minX = verts[0].x;
            maxX = verts[0].x;
            minY = verts[0].y;
            maxY = verts[0].y;
        }
        std::array<std::array<float, 2>, 2> store;
        for (unsigned int i(1); i < verts.size(); i++) {
 
                if (verts[i].x > maxX) {
                    maxX = verts[i].x;
                }
                if (verts[i].x < minX) {
                    minX = verts[i].x;
                }
                if (verts[i].y > maxY) {
                    maxY = verts[i].y;
                }
                if (verts[i].y < minY) {
                    minY = verts[i].y;
                }
 
 
        }
        store[0][0] = minX;
        store[0][1] = maxX;
        store[1][0] = minY;
        store[1][1] = maxY;
        return store;
    }
    /**\fn std::array<std::array<float, 3>,2> getExtends (const std::array<Vec3f, N>& verts);
    *  \brief get the minimum and the maximum x, y and z from an array of 3D vectors.
    *  \param the array of the vectors.
    *  \return an array containing the minimum and maximum values.
    */
    template <std::size_t N>
    static std::array<std::array<float, 2>,3> getExtends (const std::array<Vec3f, N>& verts) {
        float minX = 0;
        float maxX = 0;
        float minY = 0;
        float maxY = 0;
        float minZ = 0;
        float maxZ = 0;
        if (verts.size() > 0) {
            minX = verts[0].x;
            maxX = verts[0].x;
            minY = verts[0].y;
            maxY = verts[0].y;
            minZ = verts[0].z;
            maxZ = verts[0].z;
        }
        std::array<std::array<float, 2>, 3> store;
        for (unsigned int i(1); i < verts.size(); i++) {
 
 
                if (verts[i].x > maxX) {
                    maxX = verts[i].x;
                }
                if (verts[i].x < minX) {
                    minX = verts[i].x;
                }
                if (verts[i].y > maxY) {
                    maxY = verts[i].y;
                }
                if (verts[i].y < minY) {
                    minY = verts[i].y;
                }
                if (verts[i].z > maxZ) {
                    maxZ = verts[i].z;
                }
                if (verts[i].z < minZ) {
                    minZ = verts[i].z;
                }
 
 
        }
        store[0][0] = minX;
        store[0][1] = maxX;
        store[1][0] = minY;
        store[1][1] = maxY;
        store[2][0] = minZ;
        store[2][1] = maxZ;
        return store;
    }
    /**
     *\fn Vec3f getPosOnPathFromTime(Vec3f actualPos, std::vector<Vec3f> path, T time, float speed)
     *\brief check the position of an entity on a path (or on a curve) on a specified time.
     *If the time is positive, it checks the position in the future, otherwise it check the position on the past.
     *This function is often used for movement prediction, networking corrections and artificial intelligence.
     *\param
     * Vec3f actualPos : the actual position of the entity.
     * std::vector<Vec3f> path : the points of the path or of the curve.
     * T time : the time.
     * speed : the speed of the entity.
    */
    template <typename T>
    static Vec3f getPosOnPathFromTime(Vec3f actualPos, std::vector<Vec3f> path, T time, float speed) {
         //The number of point must be greater than 1. (Otherwise it's not a path.)
        if (path.size() > 1) {
                unsigned int currentPathIndex = 0;
                //We check the direction of the first segment of the path and the direction between out actual position and the first point of the first segment of our path.
                Vec3f dir1 = (path[currentPathIndex+1] - path[currentPathIndex]).normalize();
                Vec3f dir2 = (actualPos - path[currentPathIndex]).normalize();
                float d1 = (path[currentPathIndex+1] - path[currentPathIndex]).magnitude();
                float d2 = (actualPos - path[currentPathIndex]).magnitude();
               /* std::cout<<"dir 1 "<<dir1<<"dir 2 : "<<dir2<<"angle : "<<Math::toDegrees(dir1.getAngleBetween(dir2, dir1.cross(dir2)))<<std::endl;
                std::cout<<"dist 1 : "<<d1<<" dist 2 : "<<d2<<std::endl;*/
 
 
                //If the number of point is greater than 2, we need to check on which segment of the path our actual position is.
                //If the directions are not the same it means that the point is not on this segment of the path.
                while (currentPathIndex < path.size() - 1 && !dir2.isNulVector() && !(Math::abs(Math::toDegrees(dir1.getAngleBetween(dir2))) <= 0.5f && d2 < d1)) {
                    //We need to check if the actual position is on the next segment of the path.
                    currentPathIndex++;
                    //If we are arrived on the last point of the path, we don't need to check anymore, it means that our position is not on the path. (Or we are already at the end of the path)
                    if (currentPathIndex < path.size() - 1) {
                        dir1 = (path[currentPathIndex+1] - path[currentPathIndex]).normalize();
                        dir2 = (actualPos - path[currentPathIndex]).normalize();
                        d1 = (path[currentPathIndex+1] - path[currentPathIndex]).magnitude();
                        d2 = (actualPos - path[currentPathIndex]).magnitude();
                        /*std::cout<<"dir 1 : "<<dir1<<"dir 2 : "<<dir2<<"angle : "<<Math::toDegrees(dir1.getAngleBetween(dir2, dir1.cross(dir2)))<<std::endl;
                        std::cout<<"dist 1 : "<<d1<<" dist 2 : "<<d2<<std::endl;*/
                        //std::cout<<"current index path : "<<currentPathIndex<<std::endl;
 
                    }
                }
                //std::cout<<"current index path : "<<currentPathIndex<<"size : "<<path.size()<<std::endl;
               /*std::string str;
                std::cin>>str;*/
                //Our actual position is not on the path so we return our actual position.
                if (currentPathIndex == path.size() -1 && time > 0 || currentPathIndex == 0 && time < 0)
                    return actualPos;
                //If the time is positive, we need to check a position on the path in the future.
                if (time > 0) {
                    /*We check the direction of the segment, the next position and the distances between the two point of the segment and
                    * the distance between the next position and the first point of the segment.
                    */
                    //dir1 = (path[currentPathIndex+1] - path[currentPathIndex]).normalize();
                    Vec3f nextPos = actualPos + dir1 * speed * time;
                    dir2 = (nextPos - path[currentPathIndex]).normalize();
                    d2 = (nextPos - path[currentPathIndex]).magnitude();
                    float dot = dir1.dot(dir2);
                   // std::cout<<"actual pos : "<<actualPos<<" next pos : "<<nextPos<<std::endl;
                    //std::cout<<"dir 1 : "<<dir1<<"dir 2 : "<<dir2<<"d1 : "<<d1<<"d2 : "<<d2<<std::endl;
                    /*If the actual position is not on the same segment of the path than the final position, we need to check on which segment of the path our next position"ll be.
                    * If the distance between the two points of the segment is shorter than the distance between the next position and the first point of the path, it means that
                    * the next position is not on this segment of the path.
                    */
 
                    while (currentPathIndex < path.size() - 1 && !dir2.isNulVector() && d2 > d1) {
                        //We need to check the next position on the next segment of the path.
                        currentPathIndex++;
                        //If we are arrived on the last point of the path, we don't need to check anymore because we are at the end of the path.
                        if (currentPathIndex < path.size() - 1) {
                            dir1 = (path[currentPathIndex+1] - path[currentPathIndex]).normalize();
                            nextPos = path[currentPathIndex] + dir1 * (d2 - d1) * dot;
                            //std::cout<<"next pos : "<<nextPos<<"dir 1 : "<<dir1<<" dir 2 : "<<dir2<<" d1 : "<<d1<<" d2 : "<<d2<<"d2 - d1 : "<<d2 - d1<<" dot : "<<dot<<std::endl;
 
                            dir2 = (nextPos - path[currentPathIndex]).normalize();
                            dot = dir1.dot(dir2);
                            d1 = (path[currentPathIndex+1] - path[currentPathIndex]).magnitude();
                            d2 = (nextPos - path[currentPathIndex]).magnitude();
                            /*std::cout<<"dir 1 : "<<dir1<<" dir 2 : "<<dir2<<" d1 : "<<d1<<" d2 : "<<d2<<"d2 - d1 : "<<d2 - d1<<" dot : "<<dot<<std::endl;
                            std::string ok;
                            std::cin>>ok;*/
                        }
                    }
                    //std::cout<<"next pos : "<<nextPos<<" current path index : "<<currentPathIndex<<"size : "<<path.size()<<std::endl;
                    //If we reach the last point of the path it means that we are at the end of the path, otherwise we return the next position.
                    return (currentPathIndex == path.size() - 1)  ? path[currentPathIndex] : nextPos;
                }
                //A negative time give a position on the path in the past, so, we need to do go back on our path.
                if (time < 0) {
                    //The same as above except that we need to check the distance between the last point of the segment of the path.
                    dir1 = (path[currentPathIndex] - path[currentPathIndex+1]).normalize();
                    d1 = (path[currentPathIndex] - path[currentPathIndex+1]).magnitude();
                    Vec3f prevPos = actualPos - dir1 * speed * time;
                    dir2 = (prevPos - path[currentPathIndex+1]).normalize();
                    d2 = (prevPos - path[currentPathIndex+1]).magnitude();
                    float dot = dir1.dot(dir2);
                    while (currentPathIndex > 0 && !dir2.isNulVector() && d2 > d1) {
                        currentPathIndex--;
                        //If the current index is 0, we don't need to check anymore because we are at the beginning of the path.
                        if (currentPathIndex > 0) {
                            dir1 = (path[currentPathIndex] - path[currentPathIndex+1]).normalize();
                            prevPos = path[currentPathIndex+1] - dir1 * (d1 - d2) * dot;
                            dir2 = (prevPos - path[currentPathIndex+1]).normalize();
                            dot = dir1.dot(dir2);
                            d1 = (path[currentPathIndex+1] - path[currentPathIndex]).magnitude();
                            d2 = (prevPos - path[currentPathIndex+1]).magnitude();
                        }
                    }
                    //If we reach the last point of the path it means that we are at the end of the path, otherwise we return the next position.
                    return (currentPathIndex == 0)  ? path[currentPathIndex] : prevPos;
                }
            //If the time is equal to zero we have just to return the actual position.
            return actualPos;
        }
 
        return actualPos;
    }
    static bool overlap(Vec2f v1, Vec2f v2) {
        /*std::cout<<v1.x<<" "<<v1.y<<" "<<v2.x<<" "<<v2.y<<std::endl;
        std::cout<<((v1.y >= v2.x) || (v1.x >= v2.y))<<std::endl;*/
        return (v1.y >= v2.x) || (v1.x >= v2.y);
    }
    static Vec2f projectShapeOnAxis(Vec3f axis, std::vector<Vec3f> vertices) {
        float min, max;
        min = max = 0.f;
        if (vertices.size() > 0) {
            for (int i = 0; i < vertices.size(); i++) {
                float p = vertices[i].projOnAxis(axis);
                if (p < min) {
                    min = p;
                }
                if (p > max) {
                    max = p;
                }
            }
        }
        return Vec2f(min, max);
    }
    static Vec2f projectShapeOnAxisFromCenter(Vec3f axis, Vec3f center, std::vector<Vec3f> vertices, std::vector<Vec3f> normals) {
        Vec3f v = axis - center;
        bool found = false;
        Vec3f v1, v2, n;
        for (int i = 0; i < vertices.size() && !found; i++) {
            v1 = vertices[i] - center;
            v2 = vertices[(i+1 == vertices.size()) ? 0 : i + 1] - center;
            float angle1 = v.getAngleBetween(v1);
            float angle2 = v.getAngleBetween(v2);
            n = normals[i] - center;
            float dp = n.dot(v);
            if (dp >= 0 && (angle1 <= 0 && angle2 > 0 || angle1 > 0 && angle2 <= 0))
                found = true;
        }
        std::cout<<"v1 : "<<v1<<"v2 : "<<v2<<"n : "<<n<<"v : "<<v<<std::endl;
        float max = v.magnitude();
        float p = v.projOnAxis(n);
        float ratio = n.magnitude() / p;
        float min = max - v.magnitude() * ratio;
        std::cout<<"p : "<<p<<" n : "<<n.magnitude()<<std::endl;
        return Vec2f(min, max);
    }
    template <std::size_t N>
    static Vec2f projectShapeOnAxis(Vec3f axis, std::array<Vec3f, N> vertices) {
        float min, max;
        min = max = 0.f;
        if (vertices.size() > 0) {
            min = vertices[0].projOnAxis(axis);
            max = min;
            for (int i = 1; i < vertices.size(); i++) {
                // NOTE: the axis must be normalized to get accurate projections
                //if (!axis.isNulVector() && !vertices[i].isNulVector()) {
                    float p = vertices[i].projOnAxis(axis);
                    if (p < min) {
                        min = p;
                    }
                    if (p > max) {
                        max = p;
                    }
                //}
            }
        }
        return Vec2f(min, max);
    }
};
}
#endif

**Jedai** · 30/06/2014, 11h51

Projection de u sur n :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

(u . n) / ||n||² * n

Autrement dit, pour obtenir le projeté d'un vecteur u sur un vecteur n, on multiplie le vecteur n normalisé (divisé par sa norme) par le produit scalaire de u et n, ceci parce que si u' est le projeté de u sur n, on a :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

u . n =  u' . n = ||u'|| x ||n|| x (1 ou -1 selon que u' et n sont opposés ou non)

Autrement dit ta formule est fausse (je suppose que tu veux la longueur du projeté ?) :

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

float projection = n1.produitScalaire(v2) / n1.longueur();

Si n1 est normalisé, tu n'as même pas besoin de la division, et si tu veux vraiment la longueur, pas la distance algébrique, n'oublie pas de prendre la valeur absolue.

--
Jedaï

Invité · 30/06/2014, 21h26

Avec ta formule j'ai une projection de 16. (ce qui est encore plus loin du résultat)

Je pense finalement que cela vient d'un problème de manque de précision, je vais essayer d'utiliser des types "double" à la place des types "float".

Je crois qu'un cosinus trop imprécis peut donner des erreurs au niveau des calculs.

Peut être que j'ai un cosinus un peu trop grand ce qui me donne une projection un peu trop grande.

**bertry** · 01/07/2014, 08h09

Salut,
Si tu veux simplement déterminer si un point est dans un triangle, il y a une méthode assez simple:
http://mathworld.wolfram.com/TriangleInterior.html
Tu trouvera même une animation illustrant cette méthode ici.

Soit le triangle ABC et P, un point dont on souhaite savoir s'il est oui ou non à l'intérieur du triangle ABC:

Code :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
// Compute vectors        
v0 = C - A
v1 = B - A
v2 = P - A
 
// Compute dot products
dot00 = dot(v0, v0)
dot01 = dot(v0, v1)
dot02 = dot(v0, v2)
dot11 = dot(v1, v1)
dot12 = dot(v1, v2)
 
// Compute barycentric coordinates
invDenom = 1 / (dot00 * dot11 - dot01 * dot01)
u = (dot11 * dot02 - dot01 * dot12) * invDenom
v = (dot00 * dot12 - dot01 * dot02) * invDenom
 
// Check if point is in triangle
return (u >= 0) && (v >= 0) && (u + v < 1)

**plegat** · 01/07/2014, 13h25

Envoyé par Lolilolight

Je pense finalement que cela vient d'un problème de manque de précision, je vais essayer d'utiliser des types "double" à la place des types "float".

Je crois qu'un cosinus trop imprécis peut donner des erreurs au niveau des calculs.

Peut être que j'ai un cosinus un peu trop grand ce qui me donne une projection un peu trop grande.

A passer des double en float, on joue sur du micropouillème de degré... donc non, ton problème n'est pas là.

Sinon c'est un gros bordel de code... pas facile de trouver ton problème dans le tas (surtout que je n'ai pas trouvé ta méthode "produitScalaire" dans le colis... mais j'ai peut-être mal vu)

Serait-il possible d'avoir un code minimal représentatif stp?
Avec les coordonnées de tes vecteurs test et tout et tout...
Les formules donées par Jedai sont bonnes, donc c'est ton code qui ne l'est pas...

**Jedai** · 01/07/2014, 16h12

Envoyé par plegat

Sinon c'est un gros bordel de code... pas facile de trouver ton problème dans le tas (surtout que je n'ai pas trouvé ta méthode "produitScalaire" dans le colis... mais j'ai peut-être mal vu)

J'ai arrêté de chercher à comprendre quand j'ai vu qu'il avait trois méthodes dot(), dot2() et dot3() (dot2() est le produit scalaire).

Invité · 02/07/2014, 14h06

Re, finalement, j'ai trouvé la solution, le problème venait d'ailleurs, en fait, l'angle entre la normale et le côté ne valait pas 90°, j'ai donc rectifié cela de manière à ce que les angles entre les normales et les côtés du triangle soit toujours égaux à 90°, pour cela j'ai projeté la diagonale sur le côté et j'ai calculer le vecteur entre ce vecteur là et le centre comme ceci :

Code C++ :

Sélectionner tout - Visualiser dans une fenêtre à part

1
2
3
4
5
6
7
8
9
10
11
12
 
void BoundingPolyhedron::computeNormals() {
    normals.clear();
    for (unsigned int i = 0; i < points.size(); i++) {
        Vec3f v1 = points[i];
        Vec3f v2 = (i + 1 == points.size()) ? points[0] : points[i+1];
        Vec3f v3 = v1 - v2;
        Vec3f p = (center - v2).projOnVector(v3);
        Vec3f n = p - (center - v2);
        normals.push_back(n);
    }
}

Voilà et maintenant ça marche j'ai toujours bien un angle de 90° entre les normales et les côtés de mon triangle, du coup la projection s'effectue correctement.

Projection de vecteurs sur axes.

Mathématiques

Discussions similaires

Partager

Partager