Transférer les fichiers vers ''

tp.md

@@ -0,0 +1,116 @@
+# 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](k.png)
+
+*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](patch.png)
+
+*Visualisation d'un patch*