Code for checking xvg files in testutils.
authorDavid van der Spoel <spoel@xray.bmc.uu.se>
Mon, 27 Jul 2015 09:01:40 +0000 (11:01 +0200)
committerGerrit Code Review <gerrit@gerrit.gromacs.org>
Wed, 29 Jul 2015 19:31:53 +0000 (21:31 +0200)
Some new code to check xvg files is introduced in the testutils.
For now it only checks the contents as doubles in case they are
lines with a number of columns of numbers. The tolerance for checking
should be specified by the caller. A test case is provided.

TODO
Add checks for legends.

Change-Id: I00e3c2768c520c8bc90cf7195870a23a1ad83b76

src/gromacs/utility/stringstream.cpp
src/gromacs/utility/stringstream.h
src/testutils/CMakeLists.txt
src/testutils/cmdlinetest.cpp
src/testutils/cmdlinetest.h
src/testutils/tests/CMakeLists.txt
src/testutils/tests/xvgtest_tests.cpp [new file with mode: 0644]
src/testutils/xvgtest.cpp [new file with mode: 0644]
src/testutils/xvgtest.h [new file with mode: 0644]

index 7a93d3472c8d752fcc2c462f6d5ce04575d4d65e..77e2388ea2487cb7067790e9d13cf1e1c3b49687 100644 (file)
@@ -45,6 +45,9 @@
 
 #include <string>
 
+#include "gromacs/utility/arrayref.h"
+#include "gromacs/utility/stringutil.h"
+
 namespace gmx
 {
 
@@ -57,4 +60,40 @@ void StringOutputStream::close()
 {
 }
 
+StringInputStream::StringInputStream(const std::string &input)
+    : input_(input), pos_(0)
+{
+}
+
+StringInputStream::StringInputStream(ConstArrayRef<const char *> const &input)
+    : input_(joinStrings(input.begin(), input.end(), "\n")), pos_(0)
+{
+    input_.append("\n");
+}
+
+bool StringInputStream::readLine(std::string *line)
+{
+    if (pos_ == input_.size())
+    {
+        line->clear();
+        return false;
+    }
+    else
+    {
+        size_t newpos = input_.find("\n", pos_);
+        if (newpos == std::string::npos)
+        {
+            newpos = input_.size();
+        }
+        else
+        {
+            // To include the newline as well!
+            newpos += 1;
+        }
+        line->assign(input_.substr(pos_, newpos-pos_));
+        pos_ = newpos;
+        return true;
+    }
+}
+
 } // namespace gmx
index 899a304a12668cefe6d06de78dccdbf70e7edf62..4011688facc71121723d2768e95d41a0081ef191 100644 (file)
@@ -46,7 +46,6 @@
 
 #include <string>
 
-#include "gromacs/utility/classhelpers.h"
 #include "gromacs/utility/textstream.h"
 
 namespace gmx
@@ -75,6 +74,41 @@ class StringOutputStream : public TextOutputStream
         std::string str_;
 };
 
+template<typename T> class ConstArrayRef;
+
+/*! \libinternal \brief
+ * Helper class to convert static string data to a stream.
+ *
+ * Provides a text input stream that returns lines from a string
+ */
+class StringInputStream : public TextInputStream
+{
+    public:
+        /*! \brief
+         * Constructor that stores input lines in a string.
+         *
+         * The string is internally but no processing is done.
+         *
+         * \param[in] input String to be served by the stream.
+         */
+        explicit StringInputStream(const std::string &input);
+        /*! \brief
+         * Constructor that stores input lines in a string.
+         *
+         * The array of char * is stored as a string separated by newline.
+         *
+         * \param[in] input Array of char * to be served by the stream.
+         */
+        explicit StringInputStream(ConstArrayRef<const char *> const &input);
+
+        // From TextInputStream
+        virtual bool readLine(std::string *line);
+        virtual void close() {};
+    private:
+        std::string input_;
+        size_t      pos_;
+};
+
 } // namespace gmx
 
 #endif
index a2a5cdde513b344965bc02c57a480a67d217da30..42cf92c5c6eee2de2adbdc3df418d24f712a81fb 100644 (file)
@@ -45,7 +45,8 @@ set(TESTUTILS_SOURCES
     testfilemanager.cpp
     testfileredirector.cpp
     testinit.cpp
-    testoptions.cpp)
+    testoptions.cpp
+    xvgtest.cpp)
 
 add_library(testutils STATIC ${UNITTEST_TARGET_OPTIONS} ${TESTUTILS_SOURCES})
 set(TESTUTILS_LIBS testutils)
index 3a848c94548aeaf0ff8dd9e3e4624d20a3740334..64688b8d862d3d55a48b26d5c41330949538f215 100644 (file)
@@ -55,6 +55,7 @@
 #include "gromacs/commandline/cmdlineoptionsmodule.h"
 #include "gromacs/commandline/cmdlineprogramcontext.h"
 #include "gromacs/utility/arrayref.h"
+#include "gromacs/utility/filestream.h"
 #include "gromacs/utility/gmxassert.h"
 #include "gromacs/utility/stringutil.h"
 #include "gromacs/utility/textreader.h"
@@ -62,6 +63,7 @@
 
 #include "testutils/refdata.h"
 #include "testutils/testfilemanager.h"
+#include "testutils/xvgtest.h"
 
 namespace gmx
 {
@@ -241,10 +243,12 @@ class CommandLineTestHelper::Impl
             OutputFileInfo(const char *option, const std::string &path)
                 : option(option), path(path)
             {
+                xvg = endsWith(path, ".xvg");
             }
 
             std::string         option;
             std::string         path;
+            bool                xvg;
         };
 
         typedef std::vector<OutputFileInfo>        OutputFileList;
@@ -346,8 +350,18 @@ void CommandLineTestHelper::checkOutputFiles(TestReferenceChecker checker) const
              outfile != impl_->outputFiles_.end();
              ++outfile)
         {
-            std::string output = TextReader::readFileToString(outfile->path);
-            outputChecker.checkStringBlock(output, outfile->option.c_str());
+            if (outfile->xvg)
+            {
+                TestReferenceChecker testChecker = checker.checkCompound("File",
+                                                                         outfile->option.c_str());
+                TextInputFile        sis(outfile->path);
+                checkXvgFile(&sis, &testChecker);
+            }
+            else
+            {
+                std::string output = TextReader::readFileToString(outfile->path);
+                outputChecker.checkStringBlock(output, outfile->option.c_str());
+            }
         }
     }
 }
index 5c1feec7ea4477258598af5ac6c764790158b9a3..6d5fdc2813a1622b828ec8a5b762b114b8e03ad8 100644 (file)
@@ -287,10 +287,6 @@ class CommandLineTestHelper
          * \p filename is given to TestTemporaryFileManager to make a unique
          * filename for the temporary file, but is not otherwise used.
          *
-         * Currently, this method should not be called for an XVG file, because
-         * the comments in the beginning of the file contain timestamps and
-         * other variable information, causing the test to fail.  Best used
-         * only for custom data formats.
          */
         void setOutputFile(CommandLine *args, const char *option,
                            const char *filename);
index 351584dab7ea6fa08b6fb86b64cc3588ca669f93..6f8cb9c9b7b69a2d6725456f9623b6e3ff8963ab 100644 (file)
@@ -35,4 +35,5 @@
 gmx_add_unit_test(TestUtilsUnitTests testutils-test
                   interactivetest.cpp
                   refdata_tests.cpp
-                  testasserts_tests.cpp)
+                  testasserts_tests.cpp
+                  xvgtest_tests.cpp)
diff --git a/src/testutils/tests/xvgtest_tests.cpp b/src/testutils/tests/xvgtest_tests.cpp
new file mode 100644 (file)
index 0000000..61629ba
--- /dev/null
@@ -0,0 +1,168 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2015, by the GROMACS development team, led by
+ * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
+ * and including many others, as listed in the AUTHORS file in the
+ * top-level source directory and at http://www.gromacs.org.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+/*! \internal \file
+ * \brief
+ * Tests utilities for testing xvg files
+ *
+ * \author David van der Spoel <david.vanderspoel@icm.uu.se>
+ * \ingroup module_testutils
+ */
+#include "gmxpre.h"
+
+#include "testutils/xvgtest.h"
+
+#include <vector>
+
+#include <gtest/gtest.h>
+#include <gtest/gtest-spi.h>
+
+#include "gromacs/utility/arrayref.h"
+#include "gromacs/utility/stringstream.h"
+
+#include "testutils/refdata.h"
+#include "testutils/testasserts.h"
+
+namespace
+{
+//! Input testing data - an inline xvg file.
+const char * const input[] = {
+    "0     2905.86    -410.199",
+    "0.2     6656.67    -430.437",
+    "0.4     5262.44    -409.399",
+    "0.6     5994.69    -405.763",
+    "0.8     5941.37    -408.337",
+    "1     5869.87    -411.124"
+};
+
+TEST(XvgTests, CreateFile)
+{
+    {
+        // Create new data
+        gmx::test::TestReferenceData    data(gmx::test::erefdataUpdateAll);
+        gmx::test::TestReferenceChecker checker(data.rootChecker());
+        // Convert char array to a stream and add it to the checker
+        gmx::StringInputStream          sis(input);
+        gmx::test::checkXvgFile(&sis, &checker);
+    }
+    {
+        // Now read it back
+        gmx::test::TestReferenceData    data(gmx::test::erefdataCompare);
+        gmx::test::TestReferenceChecker checker(data.rootChecker());
+        // Convert char array to a stream and add it to the checker
+        gmx::StringInputStream          sis(input);
+        gmx::test::checkXvgFile(&sis, &checker);
+    }
+}
+
+TEST(XvgTests, CheckMissing)
+{
+    {
+        // Create new data
+        gmx::test::TestReferenceData    data(gmx::test::erefdataUpdateAll);
+        gmx::test::TestReferenceChecker checker(data.rootChecker());
+        // Convert char array to a stream and add it to the checker
+        gmx::StringInputStream          sis(input);
+        gmx::test::checkXvgFile(&sis, &checker);
+    }
+    {
+        const char * const              input[] = {
+            "0     2905.86    -410.199",
+            "0.2     6656.67    -430.437",
+            "0.4     5262.44    -409.399"
+        };
+        // Now check with missing data
+        gmx::test::TestReferenceData    data(gmx::test::erefdataCompare);
+        gmx::test::TestReferenceChecker checker(data.rootChecker());
+        gmx::StringInputStream          sis(input);
+        EXPECT_NONFATAL_FAILURE(gmx::test::checkXvgFile(&sis, &checker), "absent");
+    }
+}
+
+TEST(XvgTests, CheckExtra)
+{
+    {
+        // Create new data
+        gmx::test::TestReferenceData    data(gmx::test::erefdataUpdateAll);
+        gmx::test::TestReferenceChecker checker(data.rootChecker());
+        // Convert char array to a stream and add it to the checker
+        gmx::StringInputStream          sis(input);
+        gmx::test::checkXvgFile(&sis, &checker);
+    }
+    {
+        const char * const              input[] = {
+            "0     2905.86    -410.199",
+            "0.2     6656.67    -430.437",
+            "0.4     5262.44    -409.399",
+            "0.6     5994.69    -405.763",
+            "0.8     5941.37    -408.337",
+            "1     5869.87    -411.124",
+            "1.2     5889.87    -413.124"
+        };
+        // Now check with missing data
+        gmx::test::TestReferenceData    data(gmx::test::erefdataCompare);
+        gmx::test::TestReferenceChecker checker(data.rootChecker());
+        gmx::StringInputStream          sis(input);
+        EXPECT_NONFATAL_FAILURE(gmx::test::checkXvgFile(&sis, &checker), "Row6");
+    }
+}
+
+TEST(XvgTests, ReadIncorrect)
+{
+    {
+        // Create new data
+        gmx::test::TestReferenceData    data(gmx::test::erefdataUpdateAll);
+        gmx::test::TestReferenceChecker checker(data.rootChecker());
+        // Convert char array to a stream and add it to the checker
+        gmx::StringInputStream          sis(input);
+        gmx::test::checkXvgFile(&sis, &checker);
+    }
+    {
+        const char * const              input[] = {
+            "0     2905.86    -410.199",
+            "0.2     6656.67    -430.437",
+            "0.4     5262.44    -409.399",
+            "0.6     5994.69    -405.763",
+            "0.8     5941.37    -408.337",
+            "1     5869.87    -421.124"
+        };
+        // Now check with incorrect data
+        gmx::test::TestReferenceData    data(gmx::test::erefdataCompare);
+        gmx::test::TestReferenceChecker checker(data.rootChecker());
+        gmx::StringInputStream          sis(input);
+        EXPECT_NONFATAL_FAILURE(gmx::test::checkXvgFile(&sis, &checker), "-411");
+    }
+}
+
+} // namespace
diff --git a/src/testutils/xvgtest.cpp b/src/testutils/xvgtest.cpp
new file mode 100644 (file)
index 0000000..d13a6f2
--- /dev/null
@@ -0,0 +1,107 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2015, by the GROMACS development team, led by
+ * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
+ * and including many others, as listed in the AUTHORS file in the
+ * top-level source directory and at http://www.gromacs.org.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+/*! \internal \file
+ * \brief
+ * Implements routine to check the content of xvg files.
+ *
+ * \author David van der Spoel <david.vanderspoel@icm.uu.se>
+ * \ingroup module_testutils
+ */
+#include "gmxpre.h"
+
+#include "xvgtest.h"
+
+#include <cerrno>
+#include <cstdlib>
+
+#include <vector>
+
+#include "gromacs/utility/gmxassert.h"
+#include "gromacs/utility/stringutil.h"
+#include "gromacs/utility/textstream.h"
+
+#include "testutils/refdata.h"
+#include "testutils/testasserts.h"
+
+namespace gmx
+{
+
+namespace test
+{
+
+void checkXvgFile(TextInputStream      *input,
+                  TestReferenceChecker *checker)
+{
+    std::string line;
+    int         nrow = 0;
+
+    while (input->readLine(&line))
+    {
+        if (!((line.find("#") != line.npos) ||
+              (line.find("@") != line.npos)))
+        {
+            std::vector<std::string> split = splitString(line);
+            std::vector<real>        row;
+
+            for (std::vector<std::string>::iterator si = split.begin();
+                 (si < split.end()); ++si)
+            {
+                const char *ptr = si->c_str();
+                char       *endptr;
+                errno = 0;
+                double      dval = std::strtod(ptr, &endptr);
+                if (errno == ERANGE)
+                {
+                    GMX_THROW(InvalidInputError("Invalid value: '" + *si
+                                                + "'; it causes an overflow/underflow"));
+                }
+                if (*ptr == '\0' || *endptr != '\0')
+                {
+                    GMX_THROW(InvalidInputError("Invalid value: '" + *si
+                                                + "'; expected a number"));
+                }
+                row.push_back(dval);
+            }
+            std::string buf = formatString("Row%d", nrow++);
+            checker->checkSequence(row.begin(), row.end(), buf.c_str());
+        }
+    }
+    std::string buf = formatString("Row%d", nrow++);
+    checker->checkPresent(false, buf.c_str());
+}
+
+} // namespace test
+
+} // namespace gmx
diff --git a/src/testutils/xvgtest.h b/src/testutils/xvgtest.h
new file mode 100644 (file)
index 0000000..d41d218
--- /dev/null
@@ -0,0 +1,78 @@
+/*
+ * This file is part of the GROMACS molecular simulation package.
+ *
+ * Copyright (c) 2015, by the GROMACS development team, led by
+ * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
+ * and including many others, as listed in the AUTHORS file in the
+ * top-level source directory and at http://www.gromacs.org.
+ *
+ * GROMACS is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU Lesser General Public License
+ * as published by the Free Software Foundation; either version 2.1
+ * of the License, or (at your option) any later version.
+ *
+ * GROMACS is distributed in the hope that it will be useful,
+ * but WITHOUT ANY WARRANTY; without even the implied warranty of
+ * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+ * Lesser General Public License for more details.
+ *
+ * You should have received a copy of the GNU Lesser General Public
+ * License along with GROMACS; if not, see
+ * http://www.gnu.org/licenses, or write to the Free Software Foundation,
+ * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
+ *
+ * If you want to redistribute modifications to GROMACS, please
+ * consider that scientific software is very special. Version
+ * control is crucial - bugs must be traceable. We will be happy to
+ * consider code for inclusion in the official distribution, but
+ * derived work must not be called official GROMACS. Details are found
+ * in the README & COPYING files - if they are missing, get the
+ * official version at http://www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the research papers on the package. Check out http://www.gromacs.org.
+ */
+/*! \libinternal \file
+ * \brief
+ * Declares function to add the content of an xvg file to a checker.
+ *
+ * \author David van der Spoel <david.vanderspoel@icm.uu.se>
+ * \inlibraryapi
+ * \ingroup module_testutils
+ */
+#ifndef GMX_TESTUTILS_XVGTESTS_H
+#define GMX_TESTUTILS_XVGTESTS_H
+
+#include <string>
+
+namespace gmx
+{
+
+class TextInputStream;
+
+namespace test
+{
+
+class TestReferenceChecker;
+
+/*! \brief
+ * Adds content of xvg file to TestReferenceChecker object.
+ *
+ * A stream of strings is parsed. The columns
+ * are analyzed with a relative tolerance provided by the input.
+ * Xmgrace formatting is ignored and only multi-column data is
+ * understood.
+ *
+ * \param[in] input       Object returning the lines of the file/data
+ *                        one by one.
+ * \param[in,out] checker The checker object.
+ */
+void checkXvgFile(TextInputStream      *input,
+                  TestReferenceChecker *checker);
+
+
+} // namespace test
+
+} // namespace gmx
+
+#endif