117 lines
4.5 KiB
Markdown
117 lines
4.5 KiB
Markdown
# TP - Courbures discrètes
|
|
**Cyril Colin et Jérémy André**
|
|
|
|
## Exercice 1
|
|
Récupération des voisins d'un sommet :
|
|
Si le nombre de sommets dans le 1-voisinage est < 5, alors on ajoute les points du 2-voisinage pour calculer le patch.
|
|
```cpp
|
|
void Courbures::get_neighborhood(std::vector<MyMesh::VertexHandle> &out, const MyMesh::VertexHandle vh) {
|
|
for(...) {...} // ajout 1N à out
|
|
if(out.size() >= 5) return;
|
|
for(...) {...} // ajout 2N à out
|
|
}
|
|
```
|
|
Après avoir récupéré les voisins, la matrice de changement de base calculé pour positionner les sommets dans un repère local (pour avoir Z représentant la distance au plan XY). Puis on construit la matrice et le vecteur permettant à Eigen de calculer les coefficients du patch quadratique avec la méthode des moindres carrés. Le résultat est retourné dans une structure `QuadPatch` comprenant la matrice de changement de repère, son inverse, ainsi que les coefficients.
|
|
```cpp
|
|
QuadPatch Courbures::fit_quad(MyMesh::VertexHandle vh) {
|
|
[...]
|
|
// Récupération des voisins
|
|
get_neighborhood(neigh, vh);
|
|
[...]
|
|
|
|
[...] // Calcul de la matrice de changement de base
|
|
Eigen::Transform<double, 3, Eigen::Affine> ch_base = rot * trans;
|
|
|
|
[...]
|
|
// Calcul de la matrice / vecteur de moindres carrés linéaires (Eigen)
|
|
for(size_t i = 0; i < neigh.size(); i++) {`
|
|
MyMesh::Point point = _mesh.point(neigh[i]);
|
|
|
|
point = ch_base * point; // Application du changement de base
|
|
B[i] = -point[2]; // -zi
|
|
A(i, 0) = point[0] * point[0]; // xi^2
|
|
A(i, 1) = point[0] * point[1]; // xi*yi
|
|
A(i, 2) = point[1] * point[1]; // yi^2
|
|
A(i, 3) = point[0]; // xi
|
|
A(i, 4) = point[1]; // yi
|
|
}
|
|
|
|
// Résolution aux moindres carrés par SVD
|
|
Eigen::VectorXd coef(5);
|
|
coef = A.bdcSvd(Eigen::ComputeThinU | Eigen::ComputeThinV).solve(B);
|
|
return QuadPatch(coef, ch_base);
|
|
}
|
|
```
|
|
|
|
|
|
*Visualisation des courbures K*
|
|
|
|
## Exercice 2
|
|
Pour la seconde partie, nous avons choisi d'afficher le patch quadratique.
|
|
Avec la fonction `fit_quad` chaque sommet possèdes un attribut de type "QuadPatch" qui contient tout ce qu'il faut pour visualiser le patch.
|
|
|
|
Pour construire le maillage d'un patch sur un sommet :
|
|
- On utilise son 1-voisinage pour déterminer la taille du maillage du patch, en cherchant les coordonnées min et max.
|
|
- Puis selon un niveau de discrétisation `patch_divs`, on créé des sommets en calculant leurs positions Z avec les coefficients du patch et en les multipliant par l'inverse de la matrice de changement de base (pour retrouver sa position par rapport au maillage originel).
|
|
- Enfin, la dernière étape est d'itérer sur tous ces nouveaux sommets pour ainsi créer des triangles, permettant l'affiche du patch.
|
|
|
|
```cpp
|
|
MyMesh tesselate_quad_patch(QuadPatch q, MyMesh mesh, VertexHandle vh) {
|
|
const size_t patch_divs = 8;
|
|
|
|
MyMesh patch;
|
|
Eigen::Vector3d point(mesh.point(vh)[0], mesh.point(vh)[1], mesh.point(vh)[2]);
|
|
point = q.transform() * point;
|
|
|
|
// Recherche des dimensions englobantes du 1-anneau de vh
|
|
double xmin = point[0], xmax = point[0];
|
|
double ymin = point[1], ymax = point[1];
|
|
|
|
for (VertexHandle vi : mesh.vv_range(vh)) {
|
|
Eigen::Vector3d point(mesh.point(vi)[0], mesh.point(vi)[1], mesh.point(vi)[2]);
|
|
|
|
point = q.transform() * point;
|
|
xmin = std::min(xmin, point[0]);
|
|
xmax = std::max(xmax, point[0]);
|
|
ymin = std::min(ymin, point[1]);
|
|
ymax = std::max(ymax, point[1]);
|
|
}
|
|
|
|
// Générations des sommets
|
|
double xstep = (xmax-xmin)/static_cast<double>(patch_divs);
|
|
double ystep = (ymax-ymin)/static_cast<double>(patch_divs);
|
|
|
|
for (size_t y = 0; y < patch_divs; y++) {
|
|
for (size_t x = 0; x < patch_divs; x++) {
|
|
double dx = xmin + x * xstep;
|
|
double dy = ymin + y * ystep;
|
|
|
|
Eigen::Vector3d point(dx, dy, -q(dx, dy));
|
|
point = q.inverse_transform() * point;
|
|
patch.new_vertex(Point(point[0], point[1], point[2]));
|
|
}
|
|
}
|
|
|
|
// Générations des triangles
|
|
for (VertexHandle vhi : patch.vertices()) {
|
|
patch.set_color(vhi, MyMesh::Color(0, 1, .2));
|
|
size_t i = vhi.idx();
|
|
size_t x = i % patch_divs;
|
|
size_t y = i / patch_divs;
|
|
|
|
// On ignore la dernière colonne et dernière ligne
|
|
if (x == patch_divs - 1 || y == patch_divs - 1) continue;
|
|
|
|
patch.add_face(vhi, patch.vertex_handle(i + patch_divs), patch.vertex_handle(i + 1));
|
|
patch.add_face(patch.vertex_handle(i + 1), patch.vertex_handle(i + patch_divs),
|
|
patch.vertex_handle(i + patch_divs + 1));
|
|
}
|
|
|
|
return patch;
|
|
}
|
|
```
|
|
|
|
|
|
|
|
*Visualisation d'un patch*
|