From fa3d1e16e19ec66943347a12c920737e396a8c86 Mon Sep 17 00:00:00 2001 From: Alan Morris Date: Mon, 14 Sep 2026 16:48:24 -0600 Subject: [PATCH] Measure correspondence quality in both directions Add push (groomed to reconstruction) and an area-weighted disagreement on the groomed mesh, so tears and missing surface no longer score as perfect. Studio ranks and colors by disagreement; pull statistics are unchanged. --- Applications/shapeworks/Commands.cpp | 58 +++++- Libs/Particles/CorrespondenceEvaluation.cpp | 188 +++++++++++++++++- Libs/Particles/CorrespondenceEvaluation.h | 81 ++++++-- Libs/Python/ShapeworksPython.cpp | 20 +- Studio/Analysis/AnalysisTool.cpp | 2 +- .../Analysis/CorrespondenceQualityPanel.cpp | 116 ++++++----- Studio/Analysis/CorrespondenceQualityPanel.h | 6 +- Studio/Analysis/CorrespondenceQualityPanel.ui | 12 +- Studio/Job/CorrespondenceQualityJob.cpp | 58 ++++-- Studio/Job/CorrespondenceQualityJob.h | 24 ++- Testing/ParticlesTests/ParticlesTests.cpp | 144 +++++++++++++- docs/about/release-notes.md | 4 +- docs/studio/studio-analyze.md | 16 +- docs/workflow/analyze.md | 25 ++- 14 files changed, 610 insertions(+), 144 deletions(-) diff --git a/Applications/shapeworks/Commands.cpp b/Applications/shapeworks/Commands.cpp index f13bae816c..b05525992a 100644 --- a/Applications/shapeworks/Commands.cpp +++ b/Applications/shapeworks/Commands.cpp @@ -290,7 +290,8 @@ void CorrespondenceQualityCommand::buildParser() { const std::string desc = "Evaluate per-subject correspondence quality by reconstructing each shape from its local particles " "(biharmonic mesh warp from the cohort L1-medoid template, matching Studio's median selection) and " - "measuring distance to the groomed mesh. Reports per-subject and aggregate statistics."; + "measuring how far the reconstruction and the groomed mesh disagree, in both directions. Reports per-subject " + "and aggregate statistics."; parser.prog(prog).description(desc); parser.add_option("--name").action("store").type("string").set_default("").help( @@ -298,7 +299,8 @@ void CorrespondenceQualityCommand::buildParser() { parser.add_option("--output").action("store").type("string").set_default("").help( "Optional path to write per-subject CSV."); parser.add_option("--output_meshes").action("store").type("string").set_default("").help( - "Optional directory; write reconstructed meshes with per-vertex distance field for visual inspection."); + "Optional directory; write reconstructed meshes with the per-vertex pull distance, and groomed meshes with " + "the per-vertex disagreement and push fields, for visual inspection."); std::list methods{"point-to-cell", "point-to-point"}; parser.add_option("--method") .action("store") @@ -381,6 +383,26 @@ bool CorrespondenceQualityCommand::execute(const optparse::Values& options, Shar std::cout << " p95 = " << report.agg_norm.p95 << "\n"; std::cout << " max = " << report.agg_norm.max << "\n"; + // The pull block above keeps its original wording, which existing scripts parse, so the blocks + // below use headings that cannot be mistaken for it. + auto print_stats = [](const std::string& heading, const CorrespondenceQualityStats& raw, + const CorrespondenceQualityStats& norm) { + auto line = [](const char* label, double value, double fraction) { + std::cout << " " << label << " = " << std::setprecision(6) << value << " (" << std::setprecision(3) + << (fraction * 100.0) << "% of bbox diagonal)\n"; + }; + std::cout << heading << "\n"; + line("mean ", raw.mean, norm.mean); + line("median", raw.median, norm.median); + line("p95 ", raw.p95, norm.p95); + line("max ", raw.max, norm.max); + std::cout << std::setprecision(6); + }; + print_stats("Per-subject mean push distance (groomed -> reconstructed, area-weighted):", report.agg_push_raw, + report.agg_push_norm); + print_stats("Per-subject mean disagreement, max(push, pull) on the groomed surface (area-weighted; rank by this):", + report.agg_disagreement_raw, report.agg_disagreement_norm); + // Worst-N (sorted by normalized mean; template excluded) if (worst_n > 0) { std::vector sorted_results; @@ -390,16 +412,18 @@ bool CorrespondenceQualityCommand::execute(const optparse::Values& options, Shar } std::sort(sorted_results.begin(), sorted_results.end(), [](const CorrespondenceQualityRow& a, const CorrespondenceQualityRow& b) { - return a.norm_mean > b.norm_mean; + return a.disagreement.norm_mean > b.disagreement.norm_mean; }); const int limit = std::min(worst_n, static_cast(sorted_results.size())); - std::cout << "\nWorst " << limit << " subjects (ranked by normalized mean; template excluded):\n"; + std::cout << "\nWorst " << limit << " subjects (ranked by normalized disagreement mean; template excluded):\n"; for (int i = 0; i < limit; ++i) { const auto& r = sorted_results[i]; std::cout << " " << r.subject << " (domain " << r.domain << ")" - << " norm_mean=" << r.norm_mean << " (" << std::setprecision(3) << (r.norm_mean * 100.0) << "%)" - << std::setprecision(6) << " mean=" << r.mean_dist << " median=" << r.median_dist - << " max=" << r.max_dist << " bbox_diag=" << r.bbox_diag << "\n"; + << " disagreement=" << std::setprecision(3) << (r.disagreement.norm_mean * 100.0) << "%" + << " push=" << (r.push.norm_mean * 100.0) << "%" + << " pull=" << (r.norm_mean * 100.0) << "%" << std::setprecision(6) + << " disagreement_p99=" << r.disagreement.p99 << " disagreement_max=" << r.disagreement.max + << " bbox_diag=" << r.bbox_diag << "\n"; } } @@ -414,13 +438,29 @@ bool CorrespondenceQualityCommand::execute(const optparse::Values& options, Shar boost::filesystem::current_path(oldBasePath); return false; } + // new columns go at the end, so readers that index the original columns by position still work + auto stats_header = [](const std::string& prefix) { + std::string header; + for (const char* name : {"mean", "median", "p99", "max", "norm_mean", "norm_median", "norm_p99", "norm_max"}) { + header += "," + prefix + name; + } + return header; + }; + auto write_stats = [&csv](const CorrespondenceDistanceStats& stats) { + csv << "," << stats.mean << "," << stats.median << "," << stats.p99 << "," << stats.max << "," + << stats.norm_mean << "," << stats.norm_median << "," << stats.norm_p99 << "," << stats.norm_max; + }; csv << "subject,domain,is_template,mean_dist,median_dist,p99_dist,max_dist,bbox_diag,norm_mean,norm_median," - "norm_p99,norm_max\n"; + "norm_p99,norm_max" + << stats_header("push_") << stats_header("disagreement_") << "\n"; csv << std::fixed << std::setprecision(8); for (const auto& r : report.rows) { csv << r.subject << "," << r.domain << "," << (r.is_template ? 1 : 0) << "," << r.mean_dist << "," << r.median_dist << "," << r.p99_dist << "," << r.max_dist << "," << r.bbox_diag << "," << r.norm_mean - << "," << r.norm_median << "," << r.norm_p99 << "," << r.norm_max << "\n"; + << "," << r.norm_median << "," << r.norm_p99 << "," << r.norm_max; + write_stats(r.push); + write_stats(r.disagreement); + csv << "\n"; } SW_LOG("Wrote per-subject CSV: {}", out_path.string()); } diff --git a/Libs/Particles/CorrespondenceEvaluation.cpp b/Libs/Particles/CorrespondenceEvaluation.cpp index e4ed1a5da7..2d71bc3e82 100644 --- a/Libs/Particles/CorrespondenceEvaluation.cpp +++ b/Libs/Particles/CorrespondenceEvaluation.cpp @@ -8,6 +8,10 @@ #include #include #include +#include +#include +#include +#include #include @@ -43,6 +47,91 @@ Eigen::MatrixXd load_particles_matrix(const std::string& filename) { return m; } +//! Surface area each vertex stands for: every cell's area shared evenly among its corners. Cells are +//! fan-triangulated, so polygons other than triangles are handled too. +std::vector vertex_areas(vtkPolyData* poly_data) { + std::vector areas(poly_data->GetNumberOfPoints(), 0.0); + auto ids = vtkSmartPointer::New(); + for (vtkIdType cell = 0; cell < poly_data->GetNumberOfCells(); cell++) { + poly_data->GetCellPoints(cell, ids); + const vtkIdType corners = ids->GetNumberOfIds(); + if (corners < 3) { + continue; + } + double origin[3]; + poly_data->GetPoint(ids->GetId(0), origin); + double cell_area = 0.0; + for (vtkIdType k = 1; k + 1 < corners; k++) { + double a[3]; + double b[3]; + poly_data->GetPoint(ids->GetId(k), a); + poly_data->GetPoint(ids->GetId(k + 1), b); + double edge_a[3]; + double edge_b[3]; + vtkMath::Subtract(a, origin, edge_a); + vtkMath::Subtract(b, origin, edge_b); + double cross[3]; + vtkMath::Cross(edge_a, edge_b, cross); + cell_area += 0.5 * vtkMath::Norm(cross); + } + for (vtkIdType k = 0; k < corners; k++) { + areas[ids->GetId(k)] += cell_area / corners; + } + } + return areas; +} + +//! Mean, median, p99 and max of per-vertex values, each vertex weighted by the area it stands for, so +//! a region counts in proportion to its size however finely it happens to be meshed. +CorrespondenceDistanceStats weighted_stats(const std::vector& values, std::vector weights, + double bbox_diag) { + CorrespondenceDistanceStats stats; + if (values.empty() || weights.size() != values.size()) { + return stats; + } + + double total = std::accumulate(weights.begin(), weights.end(), 0.0); + if (!(total > 0.0)) { + // nothing to weight by (no cells), so every vertex counts the same + std::fill(weights.begin(), weights.end(), 1.0); + total = static_cast(weights.size()); + } + + std::vector order(values.size()); + std::iota(order.begin(), order.end(), 0); + std::sort(order.begin(), order.end(), [&values](size_t a, size_t b) { return values[a] < values[b]; }); + + // the smallest value with at least this fraction of the total area at or below it + auto percentile = [&](double fraction) { + const double target = fraction * total; + double cumulative = 0.0; + for (size_t i : order) { + cumulative += weights[i]; + if (cumulative >= target) { + return values[i]; + } + } + return values[order.back()]; + }; + + double weighted_sum = 0.0; + for (size_t i = 0; i < values.size(); i++) { + weighted_sum += weights[i] * values[i]; + } + stats.mean = weighted_sum / total; + stats.median = percentile(0.5); + stats.p99 = percentile(0.99); + stats.max = values[order.back()]; + + if (bbox_diag > 0.0) { + stats.norm_mean = stats.mean / bbox_diag; + stats.norm_median = stats.median / bbox_diag; + stats.norm_p99 = stats.p99 / bbox_diag; + stats.norm_max = stats.max / bbox_diag; + } + return stats; +} + } // namespace //--------------------------------------------------------------------------- @@ -61,7 +150,8 @@ CorrespondenceQualityStats CorrespondenceEvaluation::summarize(std::vector reconstructed, const Mesh& groomed, DistanceMethod method, - vtkSmartPointer* out_distance) { + vtkSmartPointer* out_distance, vtkSmartPointer* out_disagreement, + vtkSmartPointer* out_push) { CorrespondenceQualityRow row; if (!reconstructed || reconstructed->GetNumberOfPoints() == 0) { return row; @@ -70,8 +160,12 @@ CorrespondenceQualityRow CorrespondenceEvaluation::evaluate_reconstruction( const Mesh::DistanceMethod distance_method = (method == DistanceMethod::PointToPoint) ? Mesh::DistanceMethod::PointToPoint : Mesh::DistanceMethod::PointToCell; + // pull: each reconstructed vertex to the groomed surface, with the groomed cell (point-to-cell) or + // vertex (point-to-point) it is nearest Mesh recon_mesh(reconstructed); - auto field = recon_mesh.distance(groomed, distance_method)[0]; + auto pull_fields = recon_mesh.distance(groomed, distance_method); + auto field = pull_fields[0]; + auto nearest = pull_fields[1]; const int n = field->GetNumberOfTuples(); if (n == 0) { @@ -109,11 +203,68 @@ CorrespondenceQualityRow CorrespondenceEvaluation::evaluate_reconstruction( row.norm_max = row.max_dist / row.bbox_diag; } + // push: each groomed vertex to the reconstructed surface + auto push = groomed.distance(recon_mesh, distance_method)[0]; + push->SetName("push"); + + // Disagreement lives on the groomed vertices: the push distance, raised wherever a reconstructed vertex + // lands alongside with a larger pull distance. Each pull distance is carried onto the corners of the + // groomed cell that reconstructed vertex is nearest (or onto that vertex, point-to-point). There is no + // check that the two surfaces face the same way, because a fold in the reconstruction has reversed + // normals and is exactly what the pull term is there to catch. + auto groomed_poly_data = groomed.getVTKMesh(); + const vtkIdType num_groomed = groomed_poly_data->GetNumberOfPoints(); + std::vector push_values(num_groomed); + std::vector disagreement_values(num_groomed); + for (vtkIdType v = 0; v < num_groomed; v++) { + push_values[v] = std::fabs(push->GetTuple1(v)); + disagreement_values[v] = push_values[v]; + } + + auto corners = vtkSmartPointer::New(); + for (int k = 0; k < n; ++k) { + const double pull = std::fabs(field->GetTuple1(k)); + const auto target = static_cast(nearest->GetTuple1(k)); + if (distance_method == Mesh::DistanceMethod::PointToPoint) { + if (target >= 0 && target < num_groomed) { + disagreement_values[target] = std::max(disagreement_values[target], pull); + } + continue; + } + if (target < 0 || target >= groomed_poly_data->GetNumberOfCells()) { + continue; + } + groomed_poly_data->GetCellPoints(target, corners); + for (vtkIdType c = 0; c < corners->GetNumberOfIds(); c++) { + const vtkIdType corner = corners->GetId(c); + disagreement_values[corner] = std::max(disagreement_values[corner], pull); + } + } + + const auto areas = vertex_areas(groomed_poly_data); + row.push = weighted_stats(push_values, areas, row.bbox_diag); + row.disagreement = weighted_stats(disagreement_values, areas, row.bbox_diag); + if (out_distance) { field->SetName("distance"); *out_distance = field; } + if (out_disagreement) { + auto disagreement = vtkSmartPointer::New(); + disagreement->SetName("disagreement"); + disagreement->SetNumberOfComponents(1); + disagreement->SetNumberOfTuples(num_groomed); + for (vtkIdType v = 0; v < num_groomed; v++) { + disagreement->SetValue(v, disagreement_values[v]); + } + *out_disagreement = disagreement; + } + + if (out_push) { + *out_push = push; + } + return row; } @@ -121,6 +272,10 @@ CorrespondenceQualityRow CorrespondenceEvaluation::evaluate_reconstruction( void CorrespondenceEvaluation::compute_aggregates(CorrespondenceQualityReport& report) { std::vector means; std::vector norm_means; + std::vector push_means; + std::vector push_norm_means; + std::vector disagreement_means; + std::vector disagreement_norm_means; int num_template_rows = 0; for (const auto& r : report.rows) { if (r.is_template) { @@ -129,11 +284,19 @@ void CorrespondenceEvaluation::compute_aggregates(CorrespondenceQualityReport& r } means.push_back(r.mean_dist); norm_means.push_back(r.norm_mean); + push_means.push_back(r.push.mean); + push_norm_means.push_back(r.push.norm_mean); + disagreement_means.push_back(r.disagreement.mean); + disagreement_norm_means.push_back(r.disagreement.norm_mean); } report.num_template_rows = num_template_rows; report.num_evaluated = static_cast(means.size()); report.agg_raw = summarize(means); report.agg_norm = summarize(norm_means); + report.agg_push_raw = summarize(push_means); + report.agg_push_norm = summarize(push_norm_means); + report.agg_disagreement_raw = summarize(disagreement_means); + report.agg_disagreement_norm = summarize(disagreement_norm_means); } //--------------------------------------------------------------------------- @@ -265,7 +428,10 @@ CorrespondenceQualityReport CorrespondenceEvaluation::evaluate(ProjectHandle pro Mesh groomed_mesh = load_groomed_as_mesh(groomed_per_subject_domain[i][domain]); vtkSmartPointer field; - CorrespondenceQualityRow row = evaluate_reconstruction(reconstructed, groomed_mesh, method, &field); + vtkSmartPointer disagreement; + vtkSmartPointer push; + CorrespondenceQualityRow row = + evaluate_reconstruction(reconstructed, groomed_mesh, method, &field, &disagreement, &push); if (!field) continue; row.subject = name_per_subject[i]; @@ -274,14 +440,18 @@ CorrespondenceQualityReport CorrespondenceEvaluation::evaluate(ProjectHandle pro report.rows.push_back(row); if (!meshes_dir.empty()) { + const std::string stem = name_per_subject[i] + "_domain" + std::to_string(domain); + const std::string template_suffix = row.is_template ? "_TEMPLATE" : ""; + Mesh recon_mesh(reconstructed); recon_mesh.setField("distance", field, Mesh::FieldType::Point); - std::string fname = name_per_subject[i] + "_domain" + std::to_string(domain) + "_reconstructed.vtk"; - if (row.is_template) { - fname = name_per_subject[i] + "_domain" + std::to_string(domain) + "_reconstructed_TEMPLATE.vtk"; - } - boost::filesystem::path out_mesh = meshes_dir / fname; - recon_mesh.write(out_mesh.string()); + recon_mesh.write((meshes_dir / (stem + "_reconstructed" + template_suffix + ".vtk")).string()); + + // the groomed surface too, carrying the fields only it can show: a gap in the reconstruction has + // no reconstructed surface to color + groomed_mesh.setField("disagreement", disagreement, Mesh::FieldType::Point); + groomed_mesh.setField("push", push, Mesh::FieldType::Point); + groomed_mesh.write((meshes_dir / (stem + "_groomed_disagreement" + template_suffix + ".vtk")).string()); } } } diff --git a/Libs/Particles/CorrespondenceEvaluation.h b/Libs/Particles/CorrespondenceEvaluation.h index f8626c974c..a72352da1f 100644 --- a/Libs/Particles/CorrespondenceEvaluation.h +++ b/Libs/Particles/CorrespondenceEvaluation.h @@ -16,20 +16,49 @@ class Mesh; class Project; using ProjectHandle = std::shared_ptr; +//! Distance statistics over one surface, in world units and as a fraction of the subject's +//! groomed-mesh bounding box diagonal. +struct CorrespondenceDistanceStats { + double mean = 0.0; + double median = 0.0; + double p99 = 0.0; //!< the worst part of the surface, without following a single stray vertex the way max does + double max = 0.0; + double norm_mean = 0.0; + double norm_median = 0.0; + double norm_p99 = 0.0; + double norm_max = 0.0; +}; + //! Per-subject-per-domain correspondence quality result row. +/*! + * The unprefixed distances measure one direction only, from each reconstructed vertex to the groomed + * surface ("pull"). Any part of the groomed surface the reconstruction never reaches -- a tear, a + * collapsed opening, a missing appendage -- leaves every reconstructed vertex on the groomed surface, + * so pull cannot see it. `push` measures the other direction, and `disagreement` combines the two on + * the groomed surface; it is the one to rank by. + */ struct CorrespondenceQualityRow { std::string subject; int domain = 0; - double mean_dist = 0.0; //!< mean point-to-cell (or point-to-point) distance, reconstructed -> groomed - double median_dist = 0.0; //!< median per-vertex distance - double p99_dist = 0.0; //!< 99th percentile per-vertex distance, a max that ignores single outlier vertices - double max_dist = 0.0; //!< max per-vertex distance + double mean_dist = 0.0; //!< pull: mean distance from the reconstructed vertices to the groomed surface + double median_dist = 0.0; //!< pull: median per-vertex distance + double p99_dist = 0.0; //!< pull: 99th percentile per-vertex distance, a max that ignores single outlier vertices + double max_dist = 0.0; //!< pull: max per-vertex distance double bbox_diag = 0.0; //!< diagonal of the subject's groomed-mesh bounding box double norm_mean = 0.0; //!< mean_dist / bbox_diag (scale-invariant) double norm_median = 0.0; //!< median_dist / bbox_diag double norm_p99 = 0.0; //!< p99_dist / bbox_diag double norm_max = 0.0; //!< max_dist / bbox_diag bool is_template = false; //!< true for the L1-medoid template row (excluded from aggregates) + + // appended after is_template so that existing brace initializers keep their meaning + + //! from each groomed vertex to the reconstructed surface, weighted by the area each vertex stands for + CorrespondenceDistanceStats push; + + //! on each groomed vertex, the larger of its push distance and the pull distance of any reconstructed + //! vertex that lands beside it, weighted by area + CorrespondenceDistanceStats disagreement; }; //! Aggregate summary statistics. @@ -44,10 +73,14 @@ struct CorrespondenceQualityStats { struct CorrespondenceQualityReport { std::vector rows; std::string template_subject; - int num_evaluated = 0; //!< rows.size() - template rows + int num_evaluated = 0; //!< rows.size() - template rows int num_template_rows = 0; - CorrespondenceQualityStats agg_raw; //!< aggregates over raw mean_dist (template excluded) - CorrespondenceQualityStats agg_norm; //!< aggregates over bbox-normalized values (template excluded) + CorrespondenceQualityStats agg_raw; //!< aggregates over raw pull mean_dist (template excluded) + CorrespondenceQualityStats agg_norm; //!< aggregates over bbox-normalized pull (template excluded) + CorrespondenceQualityStats agg_push_raw; //!< aggregates over push.mean + CorrespondenceQualityStats agg_push_norm; //!< aggregates over push.norm_mean + CorrespondenceQualityStats agg_disagreement_raw; //!< aggregates over disagreement.mean + CorrespondenceQualityStats agg_disagreement_norm; //!< aggregates over disagreement.norm_mean }; /** @@ -56,10 +89,16 @@ struct CorrespondenceQualityReport { * * Per-subject correspondence-quality metric: reconstruct each subject's shape * from its local particles via biharmonic mesh warp from the cohort L1-medoid - * template (matches Studio's median-subject selection), then measure distance - * from the reconstruction to that subject's groomed mesh. Distances are also - * normalized by each subject's bounding-box diagonal so the metric is - * scale-invariant. + * template (matches Studio's median-subject selection), then measure how far + * that reconstruction and the subject's groomed mesh disagree. + * + * Both directions are measured. Pull, from each reconstructed vertex to the + * groomed surface, catches reconstruction that departs from the surface: folds, + * flaps, spikes. Push, from each groomed vertex to the reconstruction, catches + * groomed surface the reconstruction never reaches: tears and collapsed + * openings. Disagreement combines them into one field on the groomed mesh, the + * only surface that can show a gap. Distances are also normalized by each + * subject's bounding-box diagonal so the metric is scale-invariant. * * The template row itself is included in `rows` (with is_template=true) but * excluded from aggregate statistics — its reconstruction is near-identity @@ -79,9 +118,10 @@ class CorrespondenceEvaluation { //! (groomed, local particles) resolve. //! //! If \p output_meshes_dir is non-empty, per-subject reconstructed meshes - //! are written there as .vtk with an embedded per-vertex "distance" field. - //! The path is used verbatim (interpreted relative to the current CWD if - //! not absolute). + //! are written there as .vtk with an embedded per-vertex pull "distance" + //! field, and the groomed meshes alongside them with "disagreement" and + //! "push" fields. The path is used verbatim (interpreted relative to the + //! current CWD if not absolute). //! //! Throws std::runtime_error on setup failures (no subjects, warp failure, //! mismatched particle counts across subjects). @@ -93,18 +133,23 @@ class CorrespondenceEvaluation { //! //! Fills everything on the row except `subject`, `domain` and `is_template`, //! which the caller owns. If \p out_distance is non-null it receives the - //! per-vertex distance field (named "distance"), for surface display or - //! writing alongside the mesh. + //! pull field on the reconstruction's vertices (named "distance"). + //! \p out_disagreement and \p out_push receive fields on the groomed mesh's + //! vertices (named "disagreement" and "push"); disagreement is the field to + //! color by, since a gap in the reconstruction has no reconstructed surface + //! to show it on. //! //! Returns a default-constructed row if \p reconstructed is null or empty. static CorrespondenceQualityRow evaluate_reconstruction(vtkSmartPointer reconstructed, const Mesh& groomed, DistanceMethod method, - vtkSmartPointer* out_distance = nullptr); + vtkSmartPointer* out_distance = nullptr, + vtkSmartPointer* out_disagreement = nullptr, + vtkSmartPointer* out_push = nullptr); //! Summary statistics (mean/median/p95/max) over a set of values. static CorrespondenceQualityStats summarize(std::vector values); - //! Fill num_evaluated, num_template_rows, agg_raw and agg_norm from report.rows. + //! Fill num_evaluated, num_template_rows and the aggregates from report.rows. static void compute_aggregates(CorrespondenceQualityReport& report); }; diff --git a/Libs/Python/ShapeworksPython.cpp b/Libs/Python/ShapeworksPython.cpp index ef37807450..d477d0af6f 100644 --- a/Libs/Python/ShapeworksPython.cpp +++ b/Libs/Python/ShapeworksPython.cpp @@ -1397,6 +1397,16 @@ PYBIND11_MODULE(shapeworks_py, m) { "progress_callback"_a = nullptr, "check_abort"_a = nullptr, "surface_distance_mode"_a = false); // CorrespondenceEvaluation + py::class_(m, "CorrespondenceDistanceStats") + .def_readonly("mean", &CorrespondenceDistanceStats::mean) + .def_readonly("median", &CorrespondenceDistanceStats::median) + .def_readonly("p99", &CorrespondenceDistanceStats::p99) + .def_readonly("max", &CorrespondenceDistanceStats::max) + .def_readonly("norm_mean", &CorrespondenceDistanceStats::norm_mean) + .def_readonly("norm_median", &CorrespondenceDistanceStats::norm_median) + .def_readonly("norm_p99", &CorrespondenceDistanceStats::norm_p99) + .def_readonly("norm_max", &CorrespondenceDistanceStats::norm_max); + py::class_(m, "CorrespondenceQualityRow") .def_readonly("subject", &CorrespondenceQualityRow::subject) .def_readonly("domain", &CorrespondenceQualityRow::domain) @@ -1409,7 +1419,9 @@ PYBIND11_MODULE(shapeworks_py, m) { .def_readonly("norm_median", &CorrespondenceQualityRow::norm_median) .def_readonly("norm_p99", &CorrespondenceQualityRow::norm_p99) .def_readonly("norm_max", &CorrespondenceQualityRow::norm_max) - .def_readonly("is_template", &CorrespondenceQualityRow::is_template); + .def_readonly("is_template", &CorrespondenceQualityRow::is_template) + .def_readonly("push", &CorrespondenceQualityRow::push) + .def_readonly("disagreement", &CorrespondenceQualityRow::disagreement); py::class_(m, "CorrespondenceQualityStats") .def_readonly("mean", &CorrespondenceQualityStats::mean) @@ -1423,7 +1435,11 @@ PYBIND11_MODULE(shapeworks_py, m) { .def_readonly("num_evaluated", &CorrespondenceQualityReport::num_evaluated) .def_readonly("num_template_rows", &CorrespondenceQualityReport::num_template_rows) .def_readonly("agg_raw", &CorrespondenceQualityReport::agg_raw) - .def_readonly("agg_norm", &CorrespondenceQualityReport::agg_norm); + .def_readonly("agg_norm", &CorrespondenceQualityReport::agg_norm) + .def_readonly("agg_push_raw", &CorrespondenceQualityReport::agg_push_raw) + .def_readonly("agg_push_norm", &CorrespondenceQualityReport::agg_push_norm) + .def_readonly("agg_disagreement_raw", &CorrespondenceQualityReport::agg_disagreement_raw) + .def_readonly("agg_disagreement_norm", &CorrespondenceQualityReport::agg_disagreement_norm); py::class_ corresp_eval(m, "CorrespondenceEvaluation"); diff --git a/Studio/Analysis/AnalysisTool.cpp b/Studio/Analysis/AnalysisTool.cpp index efccfc6dfa..b93b48b1da 100644 --- a/Studio/Analysis/AnalysisTool.cpp +++ b/Studio/Analysis/AnalysisTool.cpp @@ -1769,7 +1769,7 @@ std::string AnalysisTool::get_display_feature_map() { } } - // the correspondence distance is a per-vertex field on each sample's reconstruction, so it + // the correspondence disagreement is a per-vertex field on each sample's own meshes, so it // only applies to the sample views if (correspondence_quality_panel_->get_display_distance() && (get_analysis_mode() == AnalysisTool::MODE_ALL_SAMPLES_C || diff --git a/Studio/Analysis/CorrespondenceQualityPanel.cpp b/Studio/Analysis/CorrespondenceQualityPanel.cpp index a49bbb2df4..32cfd2b13a 100644 --- a/Studio/Analysis/CorrespondenceQualityPanel.cpp +++ b/Studio/Analysis/CorrespondenceQualityPanel.cpp @@ -68,9 +68,9 @@ CorrespondenceQualityPanel::CorrespondenceQualityPanel(QWidget* parent) connect(ui_->show_distance, &QCheckBox::clicked, this, &CorrespondenceQualityPanel::show_distance_clicked); ui_->sort_metric_combo->setToolTip(StudioUtils::wrap_tooltip( - "How to rank the samples. Localized is the ratio of a sample's p99 distance to its mean: high " - "when most of the surface is fine and a small patch is badly wrong, which is what a few swapped " - "correspondence points look like.")); + "How to rank the samples by disagreement. Localized is the ratio of a sample's p99 disagreement to its " + "mean: high when most of the surface is fine and a small patch is badly wrong, which is what a few " + "swapped correspondence points look like.")); connect(ui_->sort_metric_combo, qOverload(&QComboBox::currentIndexChanged), this, &CorrespondenceQualityPanel::options_changed); connect(ui_->sort_order_combo, qOverload(&QComboBox::currentIndexChanged), this, @@ -142,19 +142,22 @@ bool CorrespondenceQualityPanel::normalized() const { return ui_->normalize_chec //--------------------------------------------------------------------------- double CorrespondenceQualityPanel::get_sort_value(const CorrespondenceQualityRow& row) const { const bool norm = normalized(); + // the panel ranks, reports and draws the two-way disagreement: the one-way pull distance cannot see + // surface the reconstruction never reaches + const auto& stats = row.disagreement; switch (ui_->sort_metric_combo->currentIndex()) { case SORT_MEDIAN: - return norm ? row.norm_median : row.median_dist; + return norm ? stats.norm_median : stats.median; case SORT_MAX: - return norm ? row.norm_max : row.max_dist; + return norm ? stats.norm_max : stats.max; case SORT_LOCALIZED: // how concentrated the error is: a few swapped particles leave most of the surface intact, // so the mean stays low while the tail spikes. p99 rather than max, which is a single // vertex and moves with one bad triangle. Scale free, so normalization does not apply. - return row.mean_dist > 0 ? row.p99_dist / row.mean_dist : 0.0; + return stats.mean > 0 ? stats.p99 / stats.mean : 0.0; case SORT_MEAN: default: - return norm ? row.norm_mean : row.mean_dist; + return norm ? stats.norm_mean : stats.mean; } } @@ -252,8 +255,8 @@ void CorrespondenceQualityPanel::show_distance_clicked() { return; } if (get_display_distance()) { - // the distance field lives on the reconstructed surfaces of each sample - session_->set_display_mode(DisplayMode::Reconstructed); + // the disagreement lives on the groomed surfaces, the only ones a gap in the reconstruction can show on + session_->set_display_mode(DisplayMode::Groomed); Q_EMIT request_samples_view(false); } session_->trigger_reinsert_shapes(); @@ -299,25 +302,31 @@ void CorrespondenceQualityPanel::handle_job_complete() { } } - // now put the measured per-vertex distances back on the surfaces - for (const auto& [shape_index, fields] : job_->get_distance_fields()) { - if (shape_index < 0 || shape_index >= static_cast(shapes.size())) { - continue; - } - auto meshes = shapes[shape_index]->get_reconstructed_meshes(true); - for (int d = 0; d < static_cast(fields.size()) && d < static_cast(meshes.meshes().size()); d++) { - auto poly_data = meshes.meshes()[d]->get_poly_data(); - if (poly_data && fields[d]) { - poly_data->GetPointData()->AddArray(fields[d]); + // now put the measured fields back on the surfaces: the pull distance on the reconstructions, and the + // disagreement on the groomed meshes, which set_point_features() leaves alone + auto apply_fields = [&shapes](const std::map>>& fields_by_shape, + DisplayMode mode) { + for (const auto& [shape_index, fields] : fields_by_shape) { + if (shape_index < 0 || shape_index >= static_cast(shapes.size())) { + continue; + } + auto meshes = shapes[shape_index]->get_meshes(mode, true); + for (int d = 0; d < static_cast(fields.size()) && d < static_cast(meshes.meshes().size()); d++) { + auto poly_data = meshes.meshes()[d]->get_poly_data(); + if (poly_data && fields[d]) { + poly_data->GetPointData()->AddArray(fields[d]); + } } } - } + }; + apply_fields(job_->get_distance_fields(), DisplayMode::Reconstructed); + apply_fields(job_->get_disagreement_fields(), DisplayMode::Groomed); ui_->show_distance->setEnabled(true); ui_->show_distance->setChecked(true); ui_->normalize_checkbox->setEnabled(true); ui_->sort_group->setEnabled(true); - session_->set_display_mode(DisplayMode::Reconstructed); + session_->set_display_mode(DisplayMode::Groomed); Q_EMIT request_samples_view(false); update_summary(); @@ -356,8 +365,8 @@ void CorrespondenceQualityPanel::update_summary() { ui_->summary_label->show(); const auto& report = job_->get_report(); - // the aggregates are always over the per-sample *mean* distance, whatever the sort metric is - const auto& stats = normalized() ? report.agg_norm : report.agg_raw; + // the aggregates are always over the per-sample *mean* disagreement, whatever the sort metric is + const auto& stats = normalized() ? report.agg_disagreement_norm : report.agg_disagreement_raw; const double scale = normalized() ? 100.0 : 1.0; const QString units = normalized() ? "% of bbox diagonal" : "world units"; @@ -373,7 +382,7 @@ void CorrespondenceQualityPanel::update_summary() { text += "Template (excluded)" + cell(template_name.toHtmlEscaped()) + ""; text += ""; - text += "

Mean distance across samples (" + units + ")

"; + text += "

Mean disagreement across samples (" + units + ")

"; text += ""; text += "" + heading("Mean") + heading("Median") + heading("p95") + heading("Max") + ""; @@ -383,8 +392,9 @@ void CorrespondenceQualityPanel::update_summary() { // these are percentiles across samples; the table's p99 column is across one sample's vertices ui_->summary_label->setToolTip(StudioUtils::wrap_tooltip( - "Distribution across samples of each sample's mean distance. The p95 here is over samples, unlike the p99 " - "column in the table, which is over the vertices of a single sample.")); + "Distribution across samples of each sample's mean disagreement, which measures the reconstruction against " + "the groomed surface in both directions. The p95 here is over samples, unlike the p99 column in the table, " + "which is over the surface of a single sample.")); ui_->summary_label->setText(text); } @@ -418,13 +428,13 @@ void CorrespondenceQualityPanel::update_table() { // the summary above reports percentiles across samples, these are across the vertices of one // sample, so say which is which rather than leaving two similar looking percentiles side by side - const QStringList tips = {"Distance from this sample's reconstruction to its groomed mesh", + const QStringList tips = {"How far this sample's reconstruction and groomed mesh disagree, in either direction", multi_domain ? "Anatomy this row measures" : QString(), - "Mean over this sample's reconstruction vertices", - "Median over this sample's reconstruction vertices", - "99th percentile of this sample's per-vertex distances: the worst part of the surface, " + "Mean over this sample's groomed surface, weighted by area", + "Median over this sample's groomed surface, weighted by area", + "99th percentile over this sample's groomed surface: the worst part of the surface, " "without following a single stray vertex the way the max does", - "Largest single per-vertex distance on this sample"}; + "Largest disagreement anywhere on this sample"}; int tip_index = 0; for (int c = 0; c < headers.size(); c++) { if (!multi_domain && tip_index == 1) { @@ -450,16 +460,25 @@ void CorrespondenceQualityPanel::update_table() { int col = 0; auto name_item = new QTableWidgetItem(name); // the column is narrow enough that most names elide, so the tooltip has to carry the full one - name_item->setToolTip(name + QString("\nbounding box diagonal: %1").arg(row.bbox_diag)); + // which direction the disagreement comes from says what kind of failure it is: pull for folds and + // flaps, push for tears and missing surface + name_item->setToolTip(name + QString("\nbounding box diagonal: %1" + "\nmean pull (reconstruction to groomed): %2%3" + "\nmean push (groomed to reconstruction): %4%3") + .arg(row.bbox_diag) + .arg((norm ? row.norm_mean : row.mean_dist) * scale) + .arg(norm ? QString(" %") : QString()) + .arg((norm ? row.push.norm_mean : row.push.mean) * scale)); table->setItem(i, col++, name_item); if (multi_domain) { table->setItem(i, col++, new QTableWidgetItem(QString::number(row.domain))); } - const double values[4] = {norm ? row.norm_mean : row.mean_dist, norm ? row.norm_median : row.median_dist, - norm ? row.norm_p99 : row.p99_dist, norm ? row.norm_max : row.max_dist}; - const double raw[4] = {row.mean_dist, row.median_dist, row.p99_dist, row.max_dist}; + const auto& stats = row.disagreement; + const double values[4] = {norm ? stats.norm_mean : stats.mean, norm ? stats.norm_median : stats.median, + norm ? stats.norm_p99 : stats.p99, norm ? stats.norm_max : stats.max}; + const double raw[4] = {stats.mean, stats.median, stats.p99, stats.max}; for (int v = 0; v < 4; v++) { auto item = new QTableWidgetItem(QString::number(values[v] * scale, 'f', 4)); item->setTextAlignment(Qt::AlignRight | Qt::AlignVCenter); // so the decimal points line up @@ -502,16 +521,17 @@ void CorrespondenceQualityPanel::update_graphs() { if (row.is_template) { // near-identity reconstruction, would flatten the rest of the chart continue; } + const auto& stats = row.disagreement; if (sorting_by_ratio()) { // the ratio is built from p99, so plot that rather than max - primary.push_back((norm ? row.norm_mean : row.mean_dist) * scale); - companion.push_back((norm ? row.norm_p99 : row.p99_dist) * scale); + primary.push_back((norm ? stats.norm_mean : stats.mean) * scale); + companion.push_back((norm ? stats.norm_p99 : stats.p99) * scale); } else if (ui_->sort_metric_combo->currentIndex() == SORT_MAX) { - primary.push_back((norm ? row.norm_max : row.max_dist) * scale); - companion.push_back((norm ? row.norm_mean : row.mean_dist) * scale); + primary.push_back((norm ? stats.norm_max : stats.max) * scale); + companion.push_back((norm ? stats.norm_mean : stats.mean) * scale); } else { primary.push_back(get_sort_value(row) * scale); - companion.push_back((norm ? row.norm_max : row.max_dist) * scale); + companion.push_back((norm ? stats.norm_max : stats.max) * scale); } } @@ -529,14 +549,14 @@ void CorrespondenceQualityPanel::update_graphs() { }; QString primary_label = ui_->sort_metric_combo->currentText(); - QString companion_label = "Max distance"; + QString companion_label = "Max disagreement"; if (sorting_by_ratio()) { - primary_label = "Mean distance"; - companion_label = "p99 distance"; + primary_label = "Mean disagreement"; + companion_label = "p99 disagreement"; } else if (sorting_by_name()) { - primary_label = "Mean distance"; + primary_label = "Mean disagreement"; } else if (ui_->sort_metric_combo->currentIndex() == SORT_MAX) { - companion_label = "Mean distance"; + companion_label = "Mean disagreement"; } std::vector series; @@ -544,10 +564,10 @@ void CorrespondenceQualityPanel::update_graphs() { series.push_back({to_vector(companion), companion_label, QColor(200, 60, 40)}); // median and p95 of the per-sample mean, matching the summary above - const auto& stats = norm ? job_->get_report().agg_norm : job_->get_report().agg_raw; - std::vector reference_lines{stats.median * scale, stats.p95 * scale}; + const auto& aggregate = norm ? job_->get_report().agg_disagreement_norm : job_->get_report().agg_disagreement_raw; + std::vector reference_lines{aggregate.median * scale, aggregate.p95 * scale}; - const QString y_label = norm ? "Distance (% of bbox diag)" : "Distance (world units)"; + const QString y_label = norm ? "Disagreement (% of bbox diag)" : "Disagreement (world units)"; const QString x_label = sorting_by_name() ? "Sample (name order)" : "Sample (table order)"; // the two series differ by more than an order of magnitude, so a linear axis would flatten the diff --git a/Studio/Analysis/CorrespondenceQualityPanel.h b/Studio/Analysis/CorrespondenceQualityPanel.h index f464daa05a..7b23eef21b 100644 --- a/Studio/Analysis/CorrespondenceQualityPanel.h +++ b/Studio/Analysis/CorrespondenceQualityPanel.h @@ -17,8 +17,8 @@ class Session; //! Panel for the correspondence quality analysis /*! * Runs CorrespondenceQualityJob and presents the result: a per-sample table that - * can be sorted by mean, median or max distance, a box plot of the distribution, - * and options to color the surface by the per-vertex distance and to reorder the + * can be sorted by mean, median or max disagreement, a box plot of the distribution, + * and options to color the groomed surface by the per-vertex disagreement and to reorder the * All Samples view worst-first so the challenging shapes come up front. */ class CorrespondenceQualityPanel : public QWidget { @@ -32,7 +32,7 @@ class CorrespondenceQualityPanel : public QWidget { void set_session(QSharedPointer session); void reset(); - //! should the per-vertex distance be shown on the surface? + //! should the per-vertex disagreement be shown on the surface? bool get_display_distance() const; //! name of the surface scalar to display diff --git a/Studio/Analysis/CorrespondenceQualityPanel.ui b/Studio/Analysis/CorrespondenceQualityPanel.ui index 3af20a4ef8..30c50bbd4d 100644 --- a/Studio/Analysis/CorrespondenceQualityPanel.ui +++ b/Studio/Analysis/CorrespondenceQualityPanel.ui @@ -214,7 +214,7 @@ QWidget#panel { - How the distance from each reconstructed vertex to the groomed surface is measured + How distances between the reconstruction and the groomed surface are measured, in both directions @@ -322,10 +322,10 @@ QWidget#panel { - Show distance on surface + Show disagreement on surface - Color each sample by its per-vertex distance to the groomed surface (reconstructed view) + Color each sample's groomed surface by where its reconstruction disagrees with it, in either direction (groomed view). The reconstructed view shows the one-way distance from the reconstruction instead. @@ -353,17 +353,17 @@ QWidget#panel { - Mean distance + Mean disagreement - Median distance + Median disagreement - Max distance + Max disagreement diff --git a/Studio/Job/CorrespondenceQualityJob.cpp b/Studio/Job/CorrespondenceQualityJob.cpp index 5eadefa2e6..54b1f01431 100644 --- a/Studio/Job/CorrespondenceQualityJob.cpp +++ b/Studio/Job/CorrespondenceQualityJob.cpp @@ -29,6 +29,7 @@ void CorrespondenceQualityJob::run() { row_shape_indices_.clear(); particle_values_.clear(); distance_fields_.clear(); + disagreement_fields_.clear(); auto shapes = session_->get_shapes(); auto non_excluded = session_->get_non_excluded_shapes(); @@ -61,6 +62,20 @@ void CorrespondenceQualityJob::run() { auto reconstructed = shape->get_reconstructed_meshes(true); auto groomed = shape->get_groomed_meshes(true); + // a previous run's fields stay on the cached meshes, so clear them before anything below can skip + // this sample, or the viewer would go on showing values measured against the old template + // copies, since MeshGroup's accessors are not const; they share the same meshes + for (auto group : {reconstructed, groomed}) { + if (!group.valid()) { + continue; + } + for (const auto& mesh : group.meshes()) { + if (mesh && mesh->get_poly_data()) { + mesh->get_poly_data()->GetPointData()->RemoveArray(FEATURE_NAME); + } + } + } + if (!reconstructed.valid() || !groomed.valid()) { SW_LOG("Correspondence quality: skipping '{}', reconstructed or groomed mesh unavailable", shape->get_display_name()); @@ -72,22 +87,31 @@ void CorrespondenceQualityJob::run() { const int num_domains = std::min(reconstructed.meshes().size(), groomed.meshes().size()); - // the distance sampled at each particle, all domains concatenated, which is the order - // Shape stores point features in + // the disagreement around each particle, all domains concatenated, which is the order Shape stores + // point features in std::vector per_particle; + // indexed by domain, with a null where a domain could not be scored, so each field lands on its own + // domain's mesh + distance_fields_[s].assign(num_domains, vtkSmartPointer()); + disagreement_fields_[s].assign(num_domains, vtkSmartPointer()); + for (int d = 0; d < num_domains; d++) { + auto particles = shape->get_particles().get_local_points(d); auto reconstructed_poly_data = reconstructed.meshes()[d]->get_poly_data(); auto groomed_poly_data = groomed.meshes()[d]->get_poly_data(); if (!reconstructed_poly_data || !groomed_poly_data) { + per_particle.insert(per_particle.end(), particles.size(), 0.0); // keep later domains aligned continue; } Mesh groomed_mesh(groomed_poly_data); vtkSmartPointer distance; + vtkSmartPointer disagreement; auto row = CorrespondenceEvaluation::evaluate_reconstruction(reconstructed_poly_data, groomed_mesh, method_, - &distance); - if (!distance) { + &distance, &disagreement); + if (!distance || !disagreement) { + per_particle.insert(per_particle.end(), particles.size(), 0.0); // keep later domains aligned continue; } @@ -97,17 +121,21 @@ void CorrespondenceQualityJob::run() { report_.rows.push_back(row); row_shape_indices_.push_back(s); - // leave the per-vertex field on the reconstructed mesh so it can be shown as a surface scalar + // leave the fields on the surfaces so they can be shown as surface scalars: the disagreement on the + // groomed mesh, the only surface a gap in the reconstruction can show on, and the one-directional + // pull distance on the reconstruction itself + disagreement->SetName(FEATURE_NAME); + groomed_poly_data->GetPointData()->AddArray(disagreement); + disagreement_fields_[s][d] = disagreement; + distance->SetName(FEATURE_NAME); reconstructed_poly_data->GetPointData()->AddArray(distance); - distance_fields_[s].push_back(distance); + distance_fields_[s][d] = distance; - // Color each particle by the error around it rather than at it. The warp inserts the - // particles into the mesh as vertices and maps them onto this shape's particles, which lie on - // its surface, so the distance at a particle is zero by construction and sampling there would - // give every glyph the same value. All the signal is in the gaps between particles, so - // assign every vertex to its nearest particle and average over that neighbourhood. - auto particles = shape->get_particles().get_local_points(d); + // Color each particle by the disagreement around it rather than at it. The reconstruction is + // warped through the particles, which lie on the groomed surface, so the two surfaces are pinned + // together at every particle and the error to look for is in the gaps between them. Assign every + // groomed vertex to its nearest particle and average over that neighbourhood. auto particle_points = vtkSmartPointer::New(); for (auto& particle : particles) { @@ -122,12 +150,12 @@ void CorrespondenceQualityJob::run() { std::vector sums(particles.size(), 0.0); std::vector counts(particles.size(), 0); - for (vtkIdType v = 0; v < reconstructed_poly_data->GetNumberOfPoints(); v++) { + for (vtkIdType v = 0; v < groomed_poly_data->GetNumberOfPoints(); v++) { double vertex[3]; - reconstructed_poly_data->GetPoint(v, vertex); + groomed_poly_data->GetPoint(v, vertex); vtkIdType id = locator->FindClosestPoint(vertex); if (id >= 0 && id < static_cast(sums.size())) { - sums[id] += std::fabs(distance->GetTuple1(v)); + sums[id] += std::fabs(disagreement->GetTuple1(v)); counts[id]++; } } diff --git a/Studio/Job/CorrespondenceQualityJob.h b/Studio/Job/CorrespondenceQualityJob.h index a6047a19a7..d52696f9fb 100644 --- a/Studio/Job/CorrespondenceQualityJob.h +++ b/Studio/Job/CorrespondenceQualityJob.h @@ -15,9 +15,10 @@ class Session; /*! * Reconstructs each sample through Studio's own configured mesh warper (the same * reconstruction shown in the viewer, using the user's chosen template and warp - * method) and measures the distance from that reconstruction back to the sample's - * groomed mesh. The per-vertex distance field is left on each reconstructed mesh - * under FEATURE_NAME so it can be displayed as a surface scalar. + * method) and measures how far that reconstruction and the sample's groomed mesh + * disagree, in both directions. The disagreement is left on each groomed mesh, and + * the one-directional pull distance on each reconstructed mesh, both under + * FEATURE_NAME so either view can display it as a surface scalar. */ class CorrespondenceQualityJob : public Job { Q_OBJECT @@ -27,7 +28,7 @@ class CorrespondenceQualityJob : public Job { void run() override; QString name() override { return "Correspondence Quality"; } - //! name of the per-vertex distance array attached to each reconstructed mesh + //! name of the per-vertex arrays: disagreement on each groomed mesh, pull distance on each reconstructed mesh static constexpr const char* FEATURE_NAME = "correspondence_distance"; const CorrespondenceQualityReport& get_report() const { return report_; } @@ -35,18 +36,24 @@ class CorrespondenceQualityJob : public Job { //! index into Session::get_shapes() for each row of the report const std::vector& get_row_shape_indices() const { return row_shape_indices_; } - //! the distance field sampled at each particle, keyed by index into Session::get_shapes(). + //! the disagreement averaged around each particle, keyed by index into Session::get_shapes(). //! The glyphs are colored from these, so they need to be applied to the shapes (on the GUI //! thread) with Shape::set_point_features() before the field can be displayed. const std::map& get_particle_values() const { return particle_values_; } - //! the per-vertex distance field for each domain of each shape, keyed by index into - //! Session::get_shapes(). Shape::set_point_features() interpolates the particle values back over - //! the mesh under the same name, so these have to be re-applied after it to survive. + //! the per-vertex pull distance on the reconstructed mesh, for each domain of each shape, keyed by index + //! into Session::get_shapes(). Shape::set_point_features() interpolates the particle values back over + //! the reconstructed mesh under the same name, so these have to be re-applied after it to survive. const std::map>>& get_distance_fields() const { return distance_fields_; } + //! the per-vertex disagreement on the groomed mesh, for each domain of each shape, keyed by index into + //! Session::get_shapes() + const std::map>>& get_disagreement_fields() const { + return disagreement_fields_; + } + private: QSharedPointer session_; CorrespondenceEvaluation::DistanceMethod method_; @@ -55,6 +62,7 @@ class CorrespondenceQualityJob : public Job { std::vector row_shape_indices_; std::map particle_values_; std::map>> distance_fields_; + std::map>> disagreement_fields_; }; } // namespace shapeworks diff --git a/Testing/ParticlesTests/ParticlesTests.cpp b/Testing/ParticlesTests/ParticlesTests.cpp index de22a222ea..8ab5d0bb4a 100644 --- a/Testing/ParticlesTests/ParticlesTests.cpp +++ b/Testing/ParticlesTests/ParticlesTests.cpp @@ -1,3 +1,5 @@ +#include + #include #include @@ -300,27 +302,30 @@ TEST(ParticlesTests, particle_normal_evaluation_test) //--------------------------------------------------------------------------- namespace { -//! flat NxN grid of triangles in the z=0 plane, spanning [0,1] in x and y -Mesh make_grid_mesh(int n) { - Eigen::MatrixXd points(n * n, 3); +//! the first `columns` columns of the NxN grid below, so the spacing is the same whatever the width +Mesh make_grid_strip_mesh(int n, int columns) { + Eigen::MatrixXd points(n * columns, 3); for (int y = 0; y < n; y++) { - for (int x = 0; x < n; x++) { - points.row(y * n + x) << static_cast(x) / (n - 1), static_cast(y) / (n - 1), 0.0; + for (int x = 0; x < columns; x++) { + points.row(y * columns + x) << static_cast(x) / (n - 1), static_cast(y) / (n - 1), 0.0; } } - Eigen::MatrixXi faces(2 * (n - 1) * (n - 1), 3); + Eigen::MatrixXi faces(2 * (n - 1) * (columns - 1), 3); int f = 0; for (int y = 0; y < n - 1; y++) { - for (int x = 0; x < n - 1; x++) { - const int i = y * n + x; - faces.row(f++) << i, i + 1, i + n; - faces.row(f++) << i + 1, i + n + 1, i + n; + for (int x = 0; x < columns - 1; x++) { + const int i = y * columns + x; + faces.row(f++) << i, i + 1, i + columns; + faces.row(f++) << i + 1, i + columns + 1, i + columns; } } return Mesh(points, faces); } +//! flat NxN grid of triangles in the z=0 plane, spanning [0,1] in x and y +Mesh make_grid_mesh(int n) { return make_grid_strip_mesh(n, n); } + } // namespace //--------------------------------------------------------------------------- @@ -411,3 +416,122 @@ TEST(CorrespondenceEvaluationTests, aggregatesExcludeTheTemplate) { ASSERT_NEAR(report.agg_raw.mean, 2.0, 1e-9); // the template row would have dominated this ASSERT_NEAR(report.agg_raw.max, 3.0, 1e-9); } + +//--------------------------------------------------------------------------- +TEST(CorrespondenceEvaluationTests, identicalMeshesHaveNoDisagreement) { + Mesh mesh = make_grid_mesh(10); + + auto row = CorrespondenceEvaluation::evaluate_reconstruction(mesh.getVTKMesh(), mesh, + CorrespondenceEvaluation::DistanceMethod::PointToCell); + + ASSERT_NEAR(row.push.max, 0.0, 1e-9); + ASSERT_NEAR(row.disagreement.max, 0.0, 1e-9); +} + +//--------------------------------------------------------------------------- +TEST(CorrespondenceEvaluationTests, uniformOffsetIsTheSameFromEitherSide) { + const double offset = 0.25; + Mesh groomed = make_grid_mesh(10); + Mesh shifted = make_grid_mesh(10); + shifted.translate(makeVector({0, 0, offset})); + + for (auto method : {CorrespondenceEvaluation::DistanceMethod::PointToCell, + CorrespondenceEvaluation::DistanceMethod::PointToPoint}) { + auto row = CorrespondenceEvaluation::evaluate_reconstruction(shifted.getVTKMesh(), groomed, method); + + ASSERT_NEAR(row.push.mean, offset, 1e-6); + ASSERT_NEAR(row.push.max, offset, 1e-6); + ASSERT_NEAR(row.disagreement.mean, offset, 1e-6); + ASSERT_NEAR(row.disagreement.norm_mean, offset / std::sqrt(2.0), 1e-6); + } +} + +//--------------------------------------------------------------------------- +// What pull alone misses: a reconstruction that covers only part of the groomed surface lies exactly +// on it, so every reconstructed vertex is at distance zero however much surface is missing. +TEST(CorrespondenceEvaluationTests, missingRegionIsCaughtByPushNotPull) { + const int n = 21; // puts a column of vertices at x = 0.5 + Mesh groomed = make_grid_mesh(n); + Mesh half = make_grid_strip_mesh(n, 11); // x <= 0.5 only + + auto row = CorrespondenceEvaluation::evaluate_reconstruction(half.getVTKMesh(), groomed, + CorrespondenceEvaluation::DistanceMethod::PointToCell); + + ASSERT_NEAR(row.max_dist, 0.0, 1e-9); // pull sees nothing wrong + + // each groomed vertex beyond x = 0.5 is (x - 0.5) from the edge of the reconstruction + ASSERT_NEAR(row.push.max, 0.5, 1e-6); + ASSERT_NEAR(row.disagreement.max, 0.5, 1e-6); + + // weighted by area over the unit square, the mean of max(0, x - 0.5) is 1/8; unweighted, the + // grid's half-area boundary vertices would pull it away from that + ASSERT_NEAR(row.push.mean, 0.125, 1e-6); + ASSERT_NEAR(row.disagreement.mean, 0.125, 1e-6); +} + +//--------------------------------------------------------------------------- +// And what push alone misses: one reconstructed vertex dragged off the surface leaves the groomed +// surface covered, so only its pull distance, carried onto the groomed vertices beneath it, sees it. +TEST(CorrespondenceEvaluationTests, offSurfaceSpikeIsCaughtByPullNotPush) { + const int n = 20; + Mesh groomed = make_grid_mesh(n); + + Mesh damaged = make_grid_mesh(n); + auto poly_data = damaged.getVTKMesh(); + double point[3]; + poly_data->GetPoint(0, point); + point[2] += 1.0; + poly_data->GetPoints()->SetPoint(0, point); + poly_data->Modified(); + + auto row = CorrespondenceEvaluation::evaluate_reconstruction(poly_data, groomed, + CorrespondenceEvaluation::DistanceMethod::PointToCell); + + ASSERT_NEAR(row.max_dist, 1.0, 1e-6); // pull sees the spike + ASSERT_LT(row.push.max, 0.1); // the groomed surface is still covered + ASSERT_NEAR(row.disagreement.max, 1.0, 1e-6); +} + +//--------------------------------------------------------------------------- +TEST(CorrespondenceEvaluationTests, disagreementFieldLivesOnTheGroomedMesh) { + const int n = 21; + Mesh groomed = make_grid_mesh(n); + Mesh half = make_grid_strip_mesh(n, 11); + + vtkSmartPointer pull; + vtkSmartPointer disagreement; + vtkSmartPointer push; + CorrespondenceEvaluation::evaluate_reconstruction( + half.getVTKMesh(), groomed, CorrespondenceEvaluation::DistanceMethod::PointToCell, &pull, &disagreement, &push); + + // pull belongs to the reconstruction's vertices; the other two to the groomed mesh's, where a gap can show + ASSERT_TRUE(pull && disagreement && push); + ASSERT_EQ(pull->GetNumberOfTuples(), static_cast(half.numPoints())); + ASSERT_EQ(disagreement->GetNumberOfTuples(), static_cast(groomed.numPoints())); + ASSERT_EQ(push->GetNumberOfTuples(), static_cast(groomed.numPoints())); + ASSERT_EQ(std::string(disagreement->GetName()), "disagreement"); + ASSERT_EQ(std::string(push->GetName()), "push"); +} + +//--------------------------------------------------------------------------- +TEST(CorrespondenceEvaluationTests, disagreementAggregatesExcludeTheTemplate) { + CorrespondenceQualityReport report; + auto add_row = [&report](const std::string& subject, double value, bool is_template) { + CorrespondenceQualityRow row; + row.subject = subject; + row.is_template = is_template; + row.push.mean = value / 2.0; + row.disagreement.mean = value; + row.disagreement.norm_mean = value / 10.0; + report.rows.push_back(row); + }; + add_row("a", 1.0, false); + add_row("b", 3.0, false); + add_row("t", 99.0, true); + + CorrespondenceEvaluation::compute_aggregates(report); + + ASSERT_NEAR(report.agg_disagreement_raw.mean, 2.0, 1e-9); + ASSERT_NEAR(report.agg_disagreement_norm.max, 0.3, 1e-9); + ASSERT_NEAR(report.agg_push_raw.mean, 1.0, 1e-9); +} diff --git a/docs/about/release-notes.md b/docs/about/release-notes.md index 26c3fd9bfd..39fbcace4d 100644 --- a/docs/about/release-notes.md +++ b/docs/about/release-notes.md @@ -14,7 +14,7 @@ * **ShapeWorks Back-end** * Registration-based particle initialization as an alternative to particle splitting: particles are spread over a single reference shape and carried onto every other shape by deformable registration (rigid → affine → SyN over distance transforms), so each shape starts optimization already holding a full set of corresponding particles (#2374) - * New `correspondence-quality` command and Python API that scores each subject by reconstructing its surface from its local particles and measuring distance back to the groomed mesh, normalized by bounding-box diagonal for comparison across anatomies. Reports mean, median, 99th percentile and max distance per subject; p99 measures the worst part of a surface without following a single stray vertex the way max does (#2612) + * New `correspondence-quality` command and Python API that scores each subject by reconstructing its surface from its local particles and measuring distance back to the groomed mesh, normalized by bounding-box diagonal for comparison across anatomies. Reports mean, median, 99th percentile and max distance per subject; p99 measures the worst part of a surface without following a single stray vertex the way max does. Distances are measured in both directions — pull, from the reconstruction to the groomed surface, and push, from the groomed surface to the reconstruction — and combined into an area-weighted disagreement on the groomed mesh, so tears and missing surface that pull alone cannot see are caught; the original pull statistics are still reported (#2612) * Large-cohort optimization speedups: the per-iteration correspondence update drops an O(P·N²) identity multiply during initialization and uses a symmetric eigensolver instead of a single-threaded general SVD, restoring multi-core use on large cohorts with no change to the result (#2574) * Geodesic remeshing fixes: `geodesic_remesh_percent` is now interpreted correctly as a percentage, geodesics are actually enabled on the remeshed surface (previously it silently fell back to Euclidean distances), and per-particle face lookups are cached against the query point (#2556) * Contours are detected by cell type rather than by inspecting the first cell, and polyline cells with more than two points are split into segments, so single-polyline contours load and optimize instead of crashing (#2377, #2457) @@ -34,7 +34,7 @@ * File → Export → Export All Meshes writes the reconstructed mesh for every subject in the project (#2281) * The glyph-size slider and auto-sizing now scale to the shape's largest dimension instead of a fixed world-unit range, so very small or very large shapes get usable glyph sizes (#2459) * Cutting-plane table edits to center and normal now take effect (#2567) - * New *Correspondence Quality* panel in the Analyze pane: scores every sample by reconstructing it from its local particles through Studio's own mesh warper — the same reconstruction shown in the viewer, using your chosen template and warp method — and measuring the distance back to that sample's groomed mesh + * New *Correspondence Quality* panel in the Analyze pane: scores every sample by reconstructing it from its local particles through Studio's own mesh warper — the same reconstruction shown in the viewer, using your chosen template and warp method — and measuring how far that reconstruction and the sample's groomed mesh disagree in both directions, then ranking by and coloring the groomed surface with that disagreement * *Show distance on surface* colors each sample's surface and particles by the per-vertex distance, so a bad region can be located and not just detected * Samples can be sorted by mean, median, p99 or max distance, or by how localized the error is, which surfaces swapped-correspondence cases whose mean distance still looks healthy; *Sort samples in view* applies the same ranking to the All Samples grid so the challenging shapes come up first, and clicking a row shows that sample on its own, or returns to all samples if it is already showing * The quality chart plots the ranked metric together with the tail of the distribution on a log axis, with median and p95 marked, so a small badly reconstructed patch stays visible even when it barely moves the mean diff --git a/docs/studio/studio-analyze.md b/docs/studio/studio-analyze.md index 213c9e836c..fc5b3c7e42 100644 --- a/docs/studio/studio-analyze.md +++ b/docs/studio/studio-analyze.md @@ -140,21 +140,23 @@ The *Particle Area Analysis* panel allows for the visualization of the area of e ## Correspondence Quality ## -The *Correspondence Quality* panel scores every sample by how well its particles describe its own surface. Each sample is reconstructed from its local particles through the mesh warper configured in the [Surface Reconstruction](surface-reconstruction.md) panel — the same reconstruction shown in the viewer, using your chosen template and warp method — and the distance from that reconstruction back to the sample's groomed mesh is measured. A sample whose particles no longer follow its own surface (a failed split, a bad initialization, an outlier shape the model does not cover) shows a large distance. +The *Correspondence Quality* panel scores every sample by how well its particles describe its own surface. Each sample is reconstructed from its local particles through the mesh warper configured in the [Surface Reconstruction](surface-reconstruction.md) panel — the same reconstruction shown in the viewer, using your chosen template and warp method — and the panel measures how far that reconstruction and the sample's groomed mesh disagree. A sample whose particles no longer follow its own surface (a failed split, a bad initialization, an outlier shape the model does not cover) shows a large disagreement. + +The disagreement is measured in both directions. From the reconstruction to the groomed surface (*pull*) catches reconstruction that departs from the surface, such as folds and flaps. From the groomed surface to the reconstruction (*push*) catches surface the reconstruction never reaches, such as a torn opening or a missing appendage, which pull alone reports as perfect because every reconstructed vertex still sits on the groomed mesh. At each groomed vertex the panel takes the larger of the two, and that disagreement is what it reports, sorts by and colors with. Note that this measures the correspondence model against the *groomed* meshes, so it reflects both optimization quality and any grooming problems upstream of it. ![ShapeWorks Studio Correspondence Quality Panel](../img/studio/studio_correspondence_quality.png) -The results are shown on the samples, so running the analysis switches you to the *Samples* tab (if you are not already on a sample view) and to the reconstructed surfaces. +The results are shown on the samples, so running the analysis switches you to the *Samples* tab (if you are not already on a sample view) and to the groomed surfaces. *Template* selects the sample everything is warped from. It is the same template the [Surface Reconstruction](surface-reconstruction.md) panel uses, so changing it here changes it there as well, and pressing *Run* rebuilds the reconstructions before measuring. *Median* picks the cohort median. Changing the template discards any results already on screen, since they were measured against the previous one. ### Reading the results -Press *Run* to compute. Distances are reported per sample as mean, median, p99 and max over the reconstruction's vertices. The p99 column is the worst part of the surface without the sensitivity to a single stray vertex that max has. With *Normalize by bounding box diagonal* checked (the default), each distance is divided by that sample's groomed bounding box diagonal and shown as a percentage, which makes samples of different size comparable. Normalization applies to the summary, the table, the chart and the sort together. +Press *Run* to compute. Disagreement is reported per sample as mean, median, p99 and max over the groomed surface, weighted by area so that a region counts in proportion to its size. The p99 column is the worst part of the surface without the sensitivity to a single stray vertex that max has. Hover over a sample's name to see its mean pull and push separately, which tells you which kind of failure it has. With *Normalize by bounding box diagonal* checked (the default), each distance is divided by that sample's groomed bounding box diagonal and shown as a percentage, which makes samples of different size comparable. Normalization applies to the summary, the table, the chart and the sort together. -The chart plots the samples in the same order as the table, on a log axis with the median and p95 marked. It draws two lines: the metric you are sorting by, and the tail of the distribution beside it — the max distance normally, p99 when sorting by *Localized* so the chart matches the ranking, and the mean when sorting by max. Two lines rather than one, because a single-line chart hides the most common failure: when a few correspondence points get swapped, only a small patch of the surface is wrong, so the mean barely moves while the tail spikes. The gap between the lines is how localized the damage is — close together means a diffusely poor reconstruction, a wide gap means a small bad region on an otherwise good one. +The chart plots the samples in the same order as the table, on a log axis with the median and p95 marked. It draws two lines: the metric you are sorting by, and the tail of the distribution beside it — the max disagreement normally, p99 when sorting by *Localized* so the chart matches the ranking, and the mean when sorting by max. Two lines rather than one, because a single-line chart hides the most common failure: when a few correspondence points get swapped, only a small patch of the surface is wrong, so the mean barely moves while the tail spikes. The gap between the lines is how localized the damage is — close together means a diffusely poor reconstruction, a wide gap means a small bad region on an otherwise good one. For the same reason *Sort by* offers **Localized (p99 / mean)**, which ranks by how concentrated each sample's error is rather than how large it is, bringing swapped-particle cases to the top even when their mean distance looks healthy. A high ratio means most of the surface is fine and a small patch is badly wrong; a low one means the error is spread evenly and the reconstruction is uniformly mediocre. p99 rather than max, so one stray vertex cannot push a sample up the ranking, and the ratio is already scale-free, so the normalize option does not change it. @@ -164,13 +166,13 @@ Click a row to show that sample on its own in the viewer, which is the quickest ### Finding the challenging shapes -Use *Sort by* to rank the table by mean, median, p99 or max distance, by how localized the error is, or by name, in descending order (worst first) or ascending (best first). Checking *Sort samples in view* applies the same ranking to the *All Samples* view, so the most challenging shapes appear first in the grid. For a multi-domain project a sample is ranked by its worst domain. Unchecking it restores the original order. +Use *Sort by* to rank the table by mean, median, p99 or max disagreement, by how localized the error is, or by name, in descending order (worst first) or ascending (best first). Checking *Sort samples in view* applies the same ranking to the *All Samples* view, so the most challenging shapes appear first in the grid. For a multi-domain project a sample is ranked by its worst domain. Unchecking it restores the original order. ### Seeing where it breaks down -*Show distance on surface* colors each sample by its per-vertex distance to the groomed surface, switching the view to the reconstructed surfaces where that field lives. This shows *where* correspondence breaks down, not just which samples are worst. +*Show disagreement on surface* colors each sample's groomed surface by its per-vertex disagreement, switching the view to the groomed surfaces where that field lives. This shows *where* correspondence breaks down, not just which samples are worst. The groomed surface is the only one that can show a gap: where the reconstruction tears or fails to reach an opening, there is no reconstructed surface to color. Switching to the reconstructed view shows the one-directional pull distance on the reconstruction itself instead, which places folds and flaps where they occur. -The particles are colored by the average distance over the part of the surface nearest to each one, rather than by the distance at the particle itself. The reconstruction is warped through the particles and the particles lie on the surface, so the distance at a particle is zero whatever the state of the model; the error to look for is always in the gaps between them. +The particles are colored by the average disagreement over the part of the groomed surface nearest to each one, rather than by the value at the particle itself. The reconstruction is warped through the particles and the particles lie on the surface, so the two surfaces are pinned together at every particle whatever the state of the model; the error to look for is always in the gaps between them. ![ShapeWorks Studio Correspondence Quality Surface Distance](../img/studio/studio_correspondence_quality_surface.png) diff --git a/docs/workflow/analyze.md b/docs/workflow/analyze.md index 25e4fcb8b4..9a203afa24 100644 --- a/docs/workflow/analyze.md +++ b/docs/workflow/analyze.md @@ -19,7 +19,15 @@ You can scroll through the dataset and zoom in and out to inspect fewer or more Scrolling through the model tells you whether correspondence *looks* right. To quantify it per subject, use the `correspondence-quality` command, or the [Correspondence Quality panel](../studio/studio-analyze.md#correspondence-quality) in Studio, which additionally lets you sort the samples worst-first so the challenging shapes come up front. -For each subject, ShapeWorks reconstructs the surface from that subject's **local** particles by biharmonic mesh warp from the cohort template, and then measures the distance from the reconstruction back to that subject's groomed mesh. A subject whose particles no longer describe its own surface — a failed split, a bad initialization, an outlier shape the model does not cover — shows up as a large distance. +For each subject, ShapeWorks reconstructs the surface from that subject's **local** particles by biharmonic mesh warp from the cohort template, and then measures how far the reconstruction and that subject's groomed mesh disagree. A subject whose particles no longer describe its own surface — a failed split, a bad initialization, an outlier shape the model does not cover — shows up as a large disagreement. + +The disagreement is measured in both directions, because each one is blind to a different failure: + +* **Pull** measures from each reconstructed vertex to the groomed surface. It catches reconstruction that departs from the surface — folds, flaps, spikes — but not surface the reconstruction never reaches: a torn opening or a missing appendage leaves every reconstructed vertex sitting on the groomed mesh, so pull reports nothing wrong. +* **Push** measures from each groomed vertex to the reconstruction. It catches those gaps, but not a flap that sticks out while the groomed surface stays covered. +* **Disagreement** combines the two on the groomed mesh: at each groomed vertex, the larger of its push distance and the pull distance of any reconstructed vertex that lands beside it. Its area-weighted mean is the number to rank by, and it is the field to color by, since the groomed mesh is the only surface on which a gap can be shown. + +The original pull statistics are still reported, unchanged, alongside push and disagreement. The template is the cohort L1-medoid (the same subject Studio picks as the median shape). Its own reconstruction is near-identity, so its row is reported but excluded from the aggregate statistics, which would otherwise be skewed on small cohorts. @@ -39,13 +47,13 @@ shapeworks correspondence-quality --name | --- | --- | | `--name` | Path to the project file (`.swproj` or `.xlsx`). Required. | | `--output` | Write the per-subject table to CSV. | -| `--output_meshes` | Write each reconstructed mesh as `.vtk` with a per-vertex `distance` field, for visual inspection of *where* correspondence breaks down. | +| `--output_meshes` | Write each reconstructed mesh as `_domain_reconstructed.vtk` with the per-vertex pull `distance` field, and each groomed mesh as `_domain_groomed_disagreement.vtk` with per-vertex `disagreement` and `push` fields, for visual inspection of *where* correspondence breaks down. | | `--method` | `point-to-cell` (default) or `point-to-point`. | | `--worst` | How many worst-ranked subjects to print. Default 5. | -The command prints a summary — mean, median, p95 and max of the per-subject mean distance, raw and normalized — followed by the worst-ranked subjects. The CSV has one row per subject per domain with the columns `subject`, `domain`, `is_template`, `mean_dist`, `median_dist`, `p99_dist`, `max_dist`, `bbox_diag`, `norm_mean`, `norm_median`, `norm_p99`, `norm_max`. `p99_dist` is the 99th percentile of the per-vertex distances: a measure of the worst part of the surface that, unlike `max_dist`, does not move with a single bad vertex. +The command prints a summary — mean, median, p95 and max of the per-subject mean distance, raw and normalized — for pull, push and disagreement, followed by the subjects with the worst disagreement. The CSV has one row per subject per domain. Its first columns are the pull statistics, as before: `subject`, `domain`, `is_template`, `mean_dist`, `median_dist`, `p99_dist`, `max_dist`, `bbox_diag`, `norm_mean`, `norm_median`, `norm_p99`, `norm_max`. They are followed by the same eight statistics for push and then for disagreement, prefixed `push_` and `disagreement_` (`push_mean`, `push_median`, `push_p99`, `push_max`, `push_norm_mean`, … `disagreement_norm_max`). Push and disagreement are weighted by surface area over the groomed mesh, so a region counts in proportion to its size however finely it is meshed. The p99 statistics are the 99th percentile of the per-vertex distances: a measure of the worst part of the surface that, unlike max, does not move with a single bad vertex. -*A reconstructed mesh written by `--output_meshes`, coloured by the `distance` field. Load it in Studio and select `distance` from the scalar dropdown to see where the reconstruction departs from the groomed surface.* +*A reconstructed mesh written by `--output_meshes`, coloured by the pull `distance` field. The groomed meshes written alongside carry the `disagreement` field, which also shows where the reconstruction fails to cover the groomed surface. Load either in Studio and select the field from the scalar dropdown.* ![Correspondence quality distance field](../img/workflow/correspondence_quality.png) Note that this measures the correspondence model against the *groomed* meshes, so it reflects both optimization quality and any grooming problems upstream of it. @@ -61,12 +69,17 @@ project.load("project.swproj") report = sw.CorrespondenceEvaluation.evaluate(project) print(report.template_subject, report.num_evaluated) -print(report.agg_norm.mean, report.agg_norm.p95) +print(report.agg_disagreement_norm.mean, report.agg_disagreement_norm.p95) for row in report.rows: - print(row.subject, row.domain, row.mean_dist, row.median_dist, row.p99_dist, row.max_dist, row.is_template) + print(row.subject, row.domain, row.is_template) + print(" disagreement", row.disagreement.mean, row.disagreement.p99, row.disagreement.max) + print(" push ", row.push.mean, row.push.p99, row.push.max) + print(" pull ", row.mean_dist, row.p99_dist, row.max_dist) ``` +`row.push` and `row.disagreement` each have `mean`, `median`, `p99` and `max`, plus `norm_` versions divided by the bounding-box diagonal. `report.agg_push_raw`/`agg_push_norm` and `report.agg_disagreement_raw`/`agg_disagreement_norm` summarize them across subjects the way `agg_raw`/`agg_norm` summarize pull. + `evaluate()` also takes `method` (`sw.CorrespondenceEvaluation.DistanceMethod.PointToCell` or `PointToPoint`) and `output_meshes_dir`. The project's relative paths are resolved against the current working directory, so run from the project's directory. ## Running ShapeWorks Studio