From 5aed3d886b679d11b373954d14489d28fcad286d Mon Sep 17 00:00:00 2001
From: Greg Landrum <greg.landrum@gmail.com>
Date: Mon, 16 Dec 2024 16:59:13 +0100
Subject: [PATCH] Fixes #8096 (#8104)

* add test

* Fixes #8096

* change in response to review
---
 External/pubchem_shape/PubChemShape.cpp |  7 ++++---
 External/pubchem_shape/test.cpp         | 15 +++++++++++++++
 2 files changed, 19 insertions(+), 3 deletions(-)

diff --git a/External/pubchem_shape/PubChemShape.cpp b/External/pubchem_shape/PubChemShape.cpp
index 5bb3a76df..014fa3676 100644
--- a/External/pubchem_shape/PubChemShape.cpp
+++ b/External/pubchem_shape/PubChemShape.cpp
@@ -263,16 +263,17 @@ ShapeInput PrepareConformer(const ROMol &mol, int confId, bool useColors) {
   res.atom_type_vector.resize(nAlignmentAtoms, 0);
 
   RDGeom::Point3D ave;
-  for (unsigned i = 0; i < nAtoms; ++i) {
+  for (unsigned i = 0, activeAtomIdx = 0; i < nAtoms; ++i) {
     unsigned int Z = mol.getAtomWithIdx(i)->getAtomicNum();
     if (Z > 1) {
       ave += conformer.getAtomPos(i);
 
-      if (vdw_radii.find(Z) == vdw_radii.end()) {
+      if (auto rad = vdw_radii.find(Z); rad != vdw_radii.end()) {
+        rad_vector[activeAtomIdx++] = rad->second;
+      } else {
         throw ValueErrorException("No VdW radius for atom with Z=" +
                                   std::to_string(Z));
       }
-      rad_vector[i] = vdw_radii.at(Z);
     }
   }
 
diff --git a/External/pubchem_shape/test.cpp b/External/pubchem_shape/test.cpp
index 6770f1df4..9dd34e290 100644
--- a/External/pubchem_shape/test.cpp
+++ b/External/pubchem_shape/test.cpp
@@ -267,3 +267,18 @@ $$$$
     CHECK_THAT(rmsd, Catch::Matchers::WithinAbs(0.017, 0.005));
   }
 }
+
+TEST_CASE("Github #8096") {
+  SECTION("as reported") {
+    auto m1 =
+        "[H]c1c([H])c([H])c([2H])c([H])c1[H] |(1.55967,1.91617,0.0546381;0.885536,1.07172,0.030849;1.38172,-0.23747,0.0274262;2.44539,-0.439501,0.0483424;0.470206,-1.27516,-0.00361916;0.856925,-2.30002,-0.00633525;-0.896665,-1.07227,-0.0310991;-1.60071,-1.87642,-0.0551085;-1.36315,0.22877,-0.0271173;-2.43593,0.379132,-0.0487835;-0.479018,1.29083,0.00359778;-0.823965,2.31421,0.00720933)|"_smiles;
+    REQUIRE(m1);
+    auto m2 =
+        "[H]c1c([H])c([H])c([H])c([H])c1[H] |(-2.06264,-0.844763,-0.0261403;-1.04035,-0.481453,-0.0114878;-0.00743655,-1.41861,-0.0137121;-0.215455,-2.47997,-0.0295909;1.29853,-0.949412,0.00507497;2.12524,-1.65277,0.00390664;1.58501,0.395878,0.0254188;2.61997,0.704365,0.0394811;0.550242,1.31385,0.0273741;0.783172,2.37039,0.0434262;-0.763786,0.88847,0.00908113;-1.60557,1.58532,0.0100194)|"_smiles;
+    REQUIRE(m2);
+    std::vector<float> matrix(12, 0.0);
+    auto [nbr_st, nbr_ct] = AlignMolecule(*m1, *m2, matrix);
+    CHECK_THAT(nbr_st, Catch::Matchers::WithinAbs(0.940, 0.005));
+    CHECK_THAT(nbr_ct, Catch::Matchers::WithinAbs(0.902, 0.005));
+  }
+}