Changeset View
Changeset View
Standalone View
Standalone View
extern/quadriflow/src/parametrizer-sing.cpp
- This file was added.
| #include "config.hpp" | |||||
| #include "field-math.hpp" | |||||
| #include "parametrizer.hpp" | |||||
| namespace qflow { | |||||
| void Parametrizer::ComputeOrientationSingularities() { | |||||
| MatrixXd &N = hierarchy.mN[0], &Q = hierarchy.mQ[0]; | |||||
| const MatrixXi& F = hierarchy.mF; | |||||
| singularities.clear(); | |||||
| for (int f = 0; f < F.cols(); ++f) { | |||||
| int index = 0; | |||||
| int abs_index = 0; | |||||
| for (int k = 0; k < 3; ++k) { | |||||
| int i = F(k, f), j = F(k == 2 ? 0 : (k + 1), f); | |||||
| auto value = | |||||
| compat_orientation_extrinsic_index_4(Q.col(i), N.col(i), Q.col(j), N.col(j)); | |||||
| index += value.second - value.first; | |||||
| abs_index += std::abs(value.second - value.first); | |||||
| } | |||||
| int index_mod = modulo(index, 4); | |||||
| if (index_mod == 1 || index_mod == 3) { | |||||
| if (index >= 4 || index < 0) { | |||||
| Q.col(F(0, f)) = -Q.col(F(0, f)); | |||||
| } | |||||
| singularities[f] = index_mod; | |||||
| } | |||||
| } | |||||
| } | |||||
| void Parametrizer::ComputePositionSingularities() { | |||||
| const MatrixXd &V = hierarchy.mV[0], &N = hierarchy.mN[0], &Q = hierarchy.mQ[0], | |||||
| &O = hierarchy.mO[0]; | |||||
| const MatrixXi& F = hierarchy.mF; | |||||
| pos_sing.clear(); | |||||
| pos_rank.resize(F.rows(), F.cols()); | |||||
| pos_index.resize(6, F.cols()); | |||||
| for (int f = 0; f < F.cols(); ++f) { | |||||
| Vector2i index = Vector2i::Zero(); | |||||
| uint32_t i0 = F(0, f), i1 = F(1, f), i2 = F(2, f); | |||||
| Vector3d q[3] = {Q.col(i0).normalized(), Q.col(i1).normalized(), Q.col(i2).normalized()}; | |||||
| Vector3d n[3] = {N.col(i0), N.col(i1), N.col(i2)}; | |||||
| Vector3d o[3] = {O.col(i0), O.col(i1), O.col(i2)}; | |||||
| Vector3d v[3] = {V.col(i0), V.col(i1), V.col(i2)}; | |||||
| int best[3]; | |||||
| double best_dp = -std::numeric_limits<double>::infinity(); | |||||
| for (int i = 0; i < 4; ++i) { | |||||
| Vector3d v0 = rotate90_by(q[0], n[0], i); | |||||
| for (int j = 0; j < 4; ++j) { | |||||
| Vector3d v1 = rotate90_by(q[1], n[1], j); | |||||
| for (int k = 0; k < 4; ++k) { | |||||
| Vector3d v2 = rotate90_by(q[2], n[2], k); | |||||
| double dp = std::min(std::min(v0.dot(v1), v1.dot(v2)), v2.dot(v0)); | |||||
| if (dp > best_dp) { | |||||
| best_dp = dp; | |||||
| best[0] = i; | |||||
| best[1] = j; | |||||
| best[2] = k; | |||||
| } | |||||
| } | |||||
| } | |||||
| } | |||||
| pos_rank(0, f) = best[0]; | |||||
| pos_rank(1, f) = best[1]; | |||||
| pos_rank(2, f) = best[2]; | |||||
| for (int k = 0; k < 3; ++k) q[k] = rotate90_by(q[k], n[k], best[k]); | |||||
| for (int k = 0; k < 3; ++k) { | |||||
| int kn = k == 2 ? 0 : (k + 1); | |||||
| double scale_x = hierarchy.mScale, scale_y = hierarchy.mScale, | |||||
| scale_x_1 = hierarchy.mScale, scale_y_1 = hierarchy.mScale; | |||||
| if (flag_adaptive_scale) { | |||||
| scale_x *= hierarchy.mS[0](0, F(k, f)); | |||||
| scale_y *= hierarchy.mS[0](1, F(k, f)); | |||||
| scale_x_1 *= hierarchy.mS[0](0, F(kn, f)); | |||||
| scale_y_1 *= hierarchy.mS[0](1, F(kn, f)); | |||||
| if (best[k] % 2 != 0) std::swap(scale_x, scale_y); | |||||
| if (best[kn] % 2 != 0) std::swap(scale_x_1, scale_y_1); | |||||
| } | |||||
| double inv_scale_x = 1.0 / scale_x, inv_scale_y = 1.0 / scale_y, | |||||
| inv_scale_x_1 = 1.0 / scale_x_1, inv_scale_y_1 = 1.0 / scale_y_1; | |||||
| std::pair<Vector2i, Vector2i> value = compat_position_extrinsic_index_4( | |||||
| v[k], n[k], q[k], o[k], v[kn], n[kn], q[kn], o[kn], scale_x, scale_y, inv_scale_x, | |||||
| inv_scale_y, scale_x_1, scale_y_1, inv_scale_x_1, inv_scale_y_1, nullptr); | |||||
| auto diff = value.first - value.second; | |||||
| index += diff; | |||||
| pos_index(k * 2, f) = diff[0]; | |||||
| pos_index(k * 2 + 1, f) = diff[1]; | |||||
| } | |||||
| if (index != Vector2i::Zero()) { | |||||
| pos_sing[f] = rshift90(index, best[0]); | |||||
| } | |||||
| } | |||||
| } | |||||
| void Parametrizer::AnalyzeValence() { | |||||
| auto& F = hierarchy.mF; | |||||
| std::map<int, int> sing; | |||||
| for (auto& f : singularities) { | |||||
| for (int i = 0; i < 3; ++i) { | |||||
| sing[F(i, f.first)] = f.second; | |||||
| } | |||||
| } | |||||
| auto& F2E = face_edgeIds; | |||||
| auto& E2E = hierarchy.mE2E; | |||||
| auto& FQ = face_edgeOrients; | |||||
| std::set<int> sing1, sing2; | |||||
| for (int i = 0; i < F2E.size(); ++i) { | |||||
| for (int j = 0; j < 3; ++j) { | |||||
| int deid = i * 3 + j; | |||||
| int sum_int = 0; | |||||
| std::vector<int> edges; | |||||
| std::vector<double> angles; | |||||
| do { | |||||
| int deid1 = deid / 3 * 3 + (deid + 2) % 3; | |||||
| deid = E2E[deid1]; | |||||
| sum_int += (FQ[deid / 3][deid % 3] + 6 - FQ[deid1 / 3][deid1 % 3]) % 4; | |||||
| } while (deid != i * 3 + j); | |||||
| if (sum_int % 4 == 2) { | |||||
| printf("OMG! valence = 2\n"); | |||||
| exit(0); | |||||
| } | |||||
| if (sum_int % 4 == 1) sing1.insert(F(j, i)); | |||||
| if (sum_int % 4 == 3) sing2.insert(F(j, i)); | |||||
| } | |||||
| } | |||||
| int count3 = 0, count4 = 0; | |||||
| for (auto& s : singularities) { | |||||
| if (s.second == 1) | |||||
| count3 += 1; | |||||
| else | |||||
| count4 += 1; | |||||
| } | |||||
| printf("singularity: <%d %d> <%d %d>\n", (int)sing1.size(), (int)sing2.size(), count3, count4); | |||||
| } | |||||
| } // namespace qflow | |||||