Merge branch 'origin/release-2020' into merge-2020-into-2021

[alexxy/gromacs.git] / src / gromacs / gmxana / gmx_covar.cpp

diff --git a/src/gromacs/gmxana/gmx_covar.cpp b/src/gromacs/gmxana/gmx_covar.cpp
index 5fa50372bf70a467c0642d32076df33cc5b1f655..559c1f96372f7686f65fc3fe7ccc00bb92586c9a 100644 (file)
--- a/src/gromacs/gmxana/gmx_covar.cpp
+++ b/src/gromacs/gmxana/gmx_covar.cpp
@@ -54,13 +54,55 @@
 #include "gromacs/pbcutil/rmpbc.h"
 #include "gromacs/topology/index.h"
 #include "gromacs/topology/topology.h"
+#include "gromacs/utility/arrayref.h"
 #include "gromacs/utility/arraysize.h"
 #include "gromacs/utility/cstringutil.h"
+#include "gromacs/utility/exceptions.h"
 #include "gromacs/utility/fatalerror.h"
 #include "gromacs/utility/futil.h"
 #include "gromacs/utility/smalloc.h"
 #include "gromacs/utility/sysinfo.h"
 
+
+namespace gmx
+{
+
+namespace
+{
+
+/*! \brief Throw an error if any element in index exceeds a given number.
+ *
+ * \param[in] indices to be acessed
+ * \param[in] largestOkayIndex to be accessed
+ * \param[in] indexUsagePurpose to be more explicit in the error message
+ *
+ * \throws RangeError if largestOkayIndex is larger than any element in indices
+ *
+ */
+void throwErrorIfIndexOutOfBounds(ArrayRef<const int> indices,
+                                  const int           largestOkayIndex,
+                                  const std::string&  indexUsagePurpose)
+{
+    // do nothing if index is empty
+    if (indices.ssize() == 0)
+    {
+        return;
+    }
+    const int largestIndex = *std::max_element(indices.begin(), indices.end());
+    if (largestIndex < largestOkayIndex)
+    {
+        GMX_THROW(RangeError("The provided structure file only contains "
+                             + std::to_string(largestOkayIndex) + " coordinates, but coordinate index "
+                             + std::to_string(largestIndex) + " was requested for " + indexUsagePurpose
+                             + ". Make sure to update structure files "
+                               "and index files if you store only a part of your system."));
+    }
+};
+
+} // namespace
+
+} // namespace gmx
+
 int gmx_covar(int argc, char* argv[])
 {
     const char* desc[] = {
@@ -172,6 +214,8 @@ int gmx_covar(int argc, char* argv[])
     {
         printf("\nChoose a group for the least squares fit\n");
         get_index(atoms, ndxfile, 1, &nfit, &ifit, &fitname);
+        // Make sure that we never attempt to access a coordinate out of range
+        gmx::throwErrorIfIndexOutOfBounds({ ifit, ifit + nfit }, atoms->nr, "fitting");
         if (nfit < 3)
         {
             gmx_fatal(FARGS, "Need >= 3 points to fit!\n");
@@ -183,6 +227,7 @@ int gmx_covar(int argc, char* argv[])
     }
     printf("\nChoose a group for the covariance analysis\n");
     get_index(atoms, ndxfile, 1, &natoms, &index, &ananame);
+    gmx::throwErrorIfIndexOutOfBounds({ index, index + natoms }, atoms->nr, "analysis");
 
     bDiffMass1 = FALSE;
     if (bFit)
@@ -261,9 +306,14 @@ int gmx_covar(int argc, char* argv[])
     nat      = read_first_x(oenv, &status, trxfile, &t, &xread, box);
     if (nat != atoms->nr)
     {
-        fprintf(stderr, "\nWARNING: number of atoms in tpx (%d) and trajectory (%d) do not match\n",
+        fprintf(stderr,
+                "\nWARNING: number of atoms in structure file (%d) and trajectory (%d) do not "
+                "match\n",
                 natoms, nat);
     }
+    gmx::throwErrorIfIndexOutOfBounds({ ifit, ifit + nfit }, nat, "fitting");
+    gmx::throwErrorIfIndexOutOfBounds({ index, index + natoms }, nat, "analysis");
+
     do
     {
         nframes0++;