Encapsulate regexp use in selections.
authorTeemu Murtola <teemu.murtola@gmail.com>
Fri, 28 Sep 2012 11:35:35 +0000 (14:35 +0300)
committerGerrit Code Review <gerrit@gerrit.gromacs.org>
Thu, 4 Oct 2012 03:29:10 +0000 (05:29 +0200)
Add a gmxregex.h header that provides a simple regular expression
matching interface, and use that in selection code instead of using
POSIX regexps directly.  Should make it easier to add regexp support
also on Windows (or at least the complexity is then within just one
file).

Also removed HAVE_SYS_TYPES_H from config.h, as sm_keywords.cpp was the
only place where it was used, and many other files included sys/types.h
unconditionally.

Change-Id: I228ad0cf200bc07a45fd745176add8ee65448789

CMakeLists.txt
src/config.h.cmakein
src/gromacs/selection/sm_keywords.cpp
src/gromacs/selection/tests/selectioncollection.cpp
src/gromacs/utility/gmxregex.cpp [new file with mode: 0644]
src/gromacs/utility/gmxregex.h [new file with mode: 0644]

index 7992721c8528fd5655f977f87f2333cf77e48fea..5001e79386fa1d8d41faa9f49d2550c20552acb9 100644 (file)
@@ -344,7 +344,6 @@ check_include_files(pwd.h        HAVE_PWD_H)
 check_include_files(pthread.h    HAVE_PTHREAD_H)
 check_include_files(dirent.h     HAVE_DIRENT_H)
 check_include_files(regex.h      HAVE_REGEX_H)
-check_include_files(sys/types.h  HAVE_SYS_TYPES_H)
 check_include_files(sys/time.h   HAVE_SYS_TIME_H)
 check_include_files(io.h                HAVE_IO_H)
 
index 8bd389a9ec6b1c576d2da06402bbb3972773ddd7..17ea8393fd1135022cbad1bf83ae88565180d9b9 100644 (file)
 /* Define to 1 if yo have the <regex.h> header file. */
 #cmakedefine HAVE_REGEX_H
 
-/* Define to 1 if you have the <sys/types.h> header file. */
-#cmakedefine HAVE_SYS_TYPES_H
-
 /* Define to 1 if you have the <sys/time.h> header file. */
 #cmakedefine HAVE_SYS_TIME_H
 
index 58fcfefa25fb90e697fd5284d7de382fa7a17cd9..9742d995a344ec6e247dea16d3b160467d5b56f5 100644 (file)
  * \author Teemu Murtola <teemu.murtola@cbr.su.se>
  * \ingroup module_selection
  */
-#ifdef HAVE_CONFIG_H
-#include <config.h>
-#endif
-
-#include <ctype.h>
-#ifdef HAVE_SYS_TYPES_H
-#include <sys/types.h>  /*old Mac needs types before regex.h*/
-#endif
-#ifdef HAVE_REGEX_H
-#include <regex.h>
-#define USE_REGEX
-#endif
+#include <cctype>
+#include <cstring>
 
 #include <string>
 
+#include <boost/shared_ptr.hpp>
+
 #include "gromacs/legacyheaders/macros.h"
 #include "gromacs/legacyheaders/smalloc.h"
 #include "gromacs/legacyheaders/string2.h"
 
 #include "gromacs/selection/selmethod.h"
 #include "gromacs/utility/exceptions.h"
+#include "gromacs/utility/gmxregex.h"
 #include "gromacs/utility/messagestringcollector.h"
 #include "gromacs/utility/stringutil.h"
 
@@ -137,35 +130,96 @@ typedef struct t_methoddata_kwreal
     real              *r;
 } t_methoddata_kwreal;
 
+namespace
+{
+
+/*! \internal \brief
+ * Single item in the list of strings/regular expressions to match.
+ *
+ * \ingroup module_selection
+ */
+class StringKeywordMatchItem
+{
+    public:
+        /*! \brief
+         * Constructs a matcher from a string.
+         *
+         * \param[in] matchType String matching type.
+         * \param[in] str       String to use for matching.
+         */
+        StringKeywordMatchItem(gmx::SelectionStringMatchType matchType,
+                               const char *str)
+            : str_(str)
+        {
+            bool bRegExp = (matchType == gmx::eStringMatchType_RegularExpression);
+            if (matchType == gmx::eStringMatchType_Auto)
+            {
+                for (size_t j = 0; j < std::strlen(str); ++j)
+                {
+                    if (std::ispunct(str[j]) && str[j] != '?' && str[j] != '*')
+                    {
+                        bRegExp = true;
+                        break;
+                    }
+                }
+            }
+            if (bRegExp)
+            {
+                if (!gmx::Regex::isSupported())
+                {
+                    GMX_THROW(gmx::InvalidInputError(gmx::formatString(
+                                    "No regular expression support, "
+                                    "cannot match \"%s\"", str)));
+                }
+                regex_.reset(new gmx::Regex(str));
+            }
+        }
+
+        /*! \brief
+         * Checks whether this item matches a string.
+         *
+         * \param[in] matchType String matching type.
+         * \param[in] value     String to match.
+         * \returns   true if this item matches \p value.
+         */
+        bool match(gmx::SelectionStringMatchType matchType,
+                   const char *value) const
+        {
+            if (matchType == gmx::eStringMatchType_Exact)
+            {
+                return str_ == value;
+            }
+            else if (regex_)
+            {
+                return gmx::regexMatch(value, *regex_);
+            }
+            else
+            {
+                return gmx_wcmatch(str_.c_str(), value) == 0;
+            }
+        }
+
+    private:
+        //! The raw string passed for the matcher.
+        std::string                     str_;
+        //! Regular expression compiled from \p str_, if applicable.
+        boost::shared_ptr<gmx::Regex>   regex_;
+};
+
 /*! \internal \brief
  * Data structure for string keyword expression evaluation.
  */
-typedef struct t_methoddata_kwstr
+struct t_methoddata_kwstr
 {
     /** Matching type for the strings. */
     gmx::SelectionStringMatchType       matchType;
     /** Array of values for the keyword. */
     char             **v;
-    /** Number of elements in the \p val array. */
-    int                n;
-    /*! \internal \brief
-     * Array of strings/regular expressions to match against.
-     */
-    struct t_methoddata_kwstr_match {
-        /** true if the expression is a regular expression, false otherwise. */
-        bool           bRegExp;
-        /** The value to match against. */
-        union {
-#ifdef USE_REGEX
-            /** Compiled regular expression if \p bRegExp is true. */
-            regex_t    r;
-#endif
-            /** The string if \p bRegExp is false; */
-            const char *s;
-        }              u;
-    }                 *m;
-    /**< Array of strings/regular expressions to match against.*/
-} t_methoddata_kwstr;
+    /** Array of strings/regular expressions to match against.*/
+    std::vector<StringKeywordMatchItem> matches;
+};
+
+} // namespace
 
 /** Parameters for integer keyword evaluation. */
 static gmx_ana_selparam_t smparams_keyword_int[] = {
@@ -450,16 +504,14 @@ evaluate_keyword_real(t_topology *top, t_trxframe *fr, t_pbc *pbc,
 /*!
  * \param[in] npar  Not used.
  * \param     param Not used.
- * \returns Pointer to the allocated data (\ref t_methoddata_kwstr).
+ * \returns Pointer to the allocated data (t_methoddata_kwstr).
  *
- * Allocates memory for a \ref t_methoddata_kwstr structure.
+ * Allocates memory for a t_methoddata_kwstr structure.
  */
 static void *
 init_data_kwstr(int npar, gmx_ana_selparam_t *param)
 {
-    t_methoddata_kwstr *data;
-
-    snew(data, 1);
+    t_methoddata_kwstr *data = new t_methoddata_kwstr();
     data->matchType = gmx::eStringMatchType_Auto;
     return data;
 }
@@ -488,83 +540,36 @@ _gmx_selelem_set_kwstr_match_type(const gmx::SelectionTreeElementPointer &sel,
  * \param[in] top   Not used.
  * \param[in] npar  Not used (should be 2).
  * \param[in] param Method parameters (should point to \ref smparams_keyword_str).
- * \param[in] data  Should point to \ref t_methoddata_kwstr.
+ * \param[in] data  Should point to t_methoddata_kwstr.
  */
 static void
 init_kwstr(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data)
 {
     t_methoddata_kwstr *d = (t_methoddata_kwstr *)data;
-    int                 i;
 
     d->v   = param[0].val.u.s;
-    d->n   = param[1].val.nr;
     /* Return if this is not the first time */
-    if (d->m)
+    if (!d->matches.empty())
     {
         return;
     }
-    snew(d->m, d->n);
-    for (i = 0; i < d->n; ++i)
+    int n = param[1].val.nr;
+    d->matches.reserve(n);
+    for (int i = 0; i < n; ++i)
     {
         const char *s = param[1].val.u.s[i];
-        bool bRegExp = (d->matchType == gmx::eStringMatchType_RegularExpression);
-        if (d->matchType == gmx::eStringMatchType_Auto)
-        {
-            for (size_t j = 0; j < strlen(s); ++j)
-            {
-                if (ispunct(s[j]) && s[j] != '?' && s[j] != '*')
-                {
-                    bRegExp = true;
-                    break;
-                }
-            }
-        }
-        if (bRegExp)
-        {
-#ifdef USE_REGEX
-            std::string buf(gmx::formatString("^%s$", s));
-            if (regcomp(&d->m[i].u.r, buf.c_str(), REG_EXTENDED | REG_NOSUB) != 0)
-            {
-                GMX_THROW(gmx::InvalidInputError(gmx::formatString(
-                                "Error in regular expression \"%s\"", s)));
-            }
-#else
-            GMX_THROW(gmx::InvalidInputError(gmx::formatString(
-                            "No regular expression support, cannot match \"%s\"", s)));
-#endif
-        }
-        else
-        {
-            d->m[i].u.s = s;
-        }
-        d->m[i].bRegExp = bRegExp;
+        d->matches.push_back(StringKeywordMatchItem(d->matchType, s));
     }
 }
 
 /*!
- * \param data Data to free (should point to a \ref t_methoddata_kwstr).
- *
- * Frees the memory allocated for t_methoddata_kwstr::val.
+ * \param data Data to free (should point to a t_methoddata_kwstr).
  */
 static void
 free_data_kwstr(void *data)
 {
     t_methoddata_kwstr *d = (t_methoddata_kwstr *)data;
-    int                 i;
-
-    for (i = 0; i < d->n; ++i)
-    {
-        if (d->m[i].bRegExp)
-        {
-#ifdef USE_REGEX
-            /* This branch should only be taken if regular expressions
-             * are available, but the ifdef is still needed. */
-            regfree(&d->m[i].u.r);
-#endif
-        }
-    }
-    sfree(d->m);
-    sfree(d);
+    delete d;
 }
 
 /*!
@@ -581,35 +586,17 @@ evaluate_keyword_str(t_topology *top, t_trxframe *fr, t_pbc *pbc,
                      gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
 {
     t_methoddata_kwstr *d = (t_methoddata_kwstr *)data;
-    int                 i, j;
-    bool                bFound;
 
     out->u.g->isize = 0;
-    for (i = 0; i < g->isize; ++i)
+    for (int i = 0; i < g->isize; ++i)
     {
-        bFound = false;
-        for (j = 0; j < d->n && !bFound; ++j)
+        for (size_t j = 0; j < d->matches.size(); ++j)
         {
-            if (d->m[j].bRegExp)
+            if (d->matches[j].match(d->matchType, d->v[i]))
             {
-#ifdef USE_REGEX
-                /* This branch should only be taken if regular expressions
-                 * are available, but the ifdef is still needed. */
-                bFound = (regexec(&d->m[j].u.r, d->v[i], 0, NULL, 0) == 0);
-#endif
+                out->u.g->index[out->u.g->isize++] = g->index[i];
+                break;
             }
-            else if (d->matchType == gmx::eStringMatchType_Exact)
-            {
-                bFound = (strcmp(d->m[j].u.s, d->v[i]) == 0);
-            }
-            else
-            {
-                bFound = (gmx_wcmatch(d->m[j].u.s, d->v[i]) == 0);
-            }
-        }
-        if (bFound)
-        {
-            out->u.g->index[out->u.g->isize++] = g->index[i];
         }
     }
 }
index 784017fb4fade7b9f4971bcd58c556cb285ab1be..e51f4ddd3fdc8dcdf6fe8b20378898c00a1e62b8 100644 (file)
  * \author Teemu Murtola <teemu.murtola@cbr.su.se>
  * \ingroup module_selection
  */
-#ifdef HAVE_CONFIG_H
-#include "config.h"
-#endif
-
 #include <gtest/gtest.h>
 
 #include "gromacs/legacyheaders/smalloc.h"
@@ -52,6 +48,7 @@
 #include "gromacs/selection/selection.h"
 #include "gromacs/utility/exceptions.h"
 #include "gromacs/utility/flags.h"
+#include "gromacs/utility/gmxregex.h"
 #include "gromacs/utility/stringutil.h"
 
 #include "testutils/refdata.h"
@@ -376,7 +373,6 @@ TEST_F(SelectionCollectionTest, ParsesSelectionsFromFile)
     EXPECT_STREQ("resname RB RC", sel_[1].selectionText());
 }
 
-#ifdef HAVE_REGEX_H
 TEST_F(SelectionCollectionTest, HandlesInvalidRegularExpressions)
 {
     ASSERT_NO_FATAL_FAILURE(loadTopology("simple.gro"));
@@ -385,16 +381,18 @@ TEST_F(SelectionCollectionTest, HandlesInvalidRegularExpressions)
             sc_.compile();
         }, gmx::InvalidInputError);
 }
-#else
+
 TEST_F(SelectionCollectionTest, HandlesUnsupportedRegularExpressions)
 {
-    ASSERT_NO_FATAL_FAILURE(loadTopology("simple.gro"));
-    EXPECT_THROW({
-            sc_.parseFromString("resname \"R[AD]\"");
-            sc_.compile();
-        }, gmx::InvalidInputError);
+    if (!gmx::Regex::isSupported())
+    {
+        ASSERT_NO_FATAL_FAILURE(loadTopology("simple.gro"));
+        EXPECT_THROW({
+                sc_.parseFromString("resname \"R[AD]\"");
+                sc_.compile();
+            }, gmx::InvalidInputError);
+    }
 }
-#endif
 
 TEST_F(SelectionCollectionTest, HandlesMissingMethodParamValue)
 {
@@ -785,7 +783,6 @@ TEST_F(SelectionCollectionDataTest, HandlesWildcardMatching)
 }
 
 
-#ifdef HAVE_REGEX_H
 TEST_F(SelectionCollectionDataTest, HandlesRegexMatching)
 {
     static const char * const selections[] = {
@@ -793,9 +790,11 @@ TEST_F(SelectionCollectionDataTest, HandlesRegexMatching)
         "resname ~ \"R[BD]\"",
         NULL
     };
-    runTest("simple.gro", selections);
+    if (gmx::Regex::isSupported())
+    {
+        runTest("simple.gro", selections);
+    }
 }
-#endif
 
 
 TEST_F(SelectionCollectionDataTest, HandlesBasicBoolean)
diff --git a/src/gromacs/utility/gmxregex.cpp b/src/gromacs/utility/gmxregex.cpp
new file mode 100644 (file)
index 0000000..2bd6cfc
--- /dev/null
@@ -0,0 +1,168 @@
+/*
+ *
+ *                This source code is part of
+ *
+ *                 G   R   O   M   A   C   S
+ *
+ *          GROningen MAchine for Chemical Simulations
+ *
+ * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
+ * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
+ * Copyright (c) 2001-2009, The GROMACS development team,
+ * check out http://www.gromacs.org for more information.
+ *
+ * This program is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU General Public License
+ * as published by the Free Software Foundation; either version 2
+ * of the License, or (at your option) any later version.
+ *
+ * If you want to redistribute modifications, 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 www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the papers on the package - you can find them in the top README file.
+ *
+ * For more info, check our website at http://www.gromacs.org
+ */
+/*! \internal \file
+ * \brief
+ * Implements regular expression wrappers.
+ *
+ * \author Teemu Murtola <teemu.murtola@cbr.su.se>
+ * \ingroup module_utility
+ */
+#include "gmxregex.h"
+
+#ifdef HAVE_CONFIG_H
+#include "config.h"
+#endif
+
+#if defined(HAVE_REGEX_H)
+// old Mac needs sys/types.h before regex.h
+#include <sys/types.h>
+#include <regex.h>
+#define USE_POSIX_REGEX
+#endif
+
+#include "exceptions.h"
+#include "stringutil.h"
+
+namespace gmx
+{
+
+// static
+bool Regex::isSupported()
+{
+#if defined(USE_POSIX_REGEX)
+    return true;
+#else
+    return false;
+#endif
+}
+
+#if defined(USE_POSIX_REGEX)
+class Regex::Impl
+{
+    public:
+        explicit Impl(const char *value)
+        {
+            compile(value);
+        }
+        explicit Impl(const std::string &value)
+        {
+            compile(value.c_str());
+        }
+        ~Impl()
+        {
+            regfree(&regex_);
+        }
+
+        bool match(const char *value) const
+        {
+            int rc = regexec(&regex_, value, 0, NULL, 0);
+            if (rc != 0 && rc != REG_NOMATCH)
+            {
+                // TODO: Handle errors.
+            }
+            return (rc == 0);
+        }
+
+    private:
+        void compile(const char *value)
+        {
+            std::string buf(formatString("^%s$", value));
+            int rc = regcomp(&regex_, buf.c_str(), REG_EXTENDED | REG_NOSUB);
+            if (rc != 0)
+            {
+                // TODO: Better error messages.
+                GMX_THROW(InvalidInputError(formatString(
+                                "Error in regular expression \"%s\"", value)));
+            }
+        }
+
+        regex_t                 regex_;
+};
+#else
+class Regex::Impl
+{
+    public:
+        explicit Impl(const char * /*value*/)
+        {
+            GMX_THROW(NotImplementedError(
+                        "Gromacs is compiled without regular expression support"));
+        }
+        explicit Impl(const std::string & /*value*/)
+        {
+            GMX_THROW(NotImplementedError(
+                        "Gromacs is compiled without regular expression support"));
+        }
+
+        bool match(const char * /*value*/) const
+        {
+            // Should never be reached.
+            GMX_THROW(NotImplementedError(
+                        "Gromacs is compiled without regular expression support"));
+        }
+};
+#endif
+
+Regex::Regex()
+    : impl_(NULL)
+{
+}
+
+Regex::Regex(const char *value)
+    : impl_(new Impl(value))
+{
+}
+
+Regex::Regex(const std::string &value)
+    : impl_(new Impl(value))
+{
+}
+
+Regex::~Regex()
+{
+}
+
+/*! \cond libapi */
+bool regexMatch(const char *str, const Regex &regex)
+{
+    if (regex.impl_.get() == NULL)
+    {
+        return false;
+    }
+    return regex.impl_->match(str);
+}
+
+bool regexMatch(const std::string &str, const Regex &regex)
+{
+    return regexMatch(str.c_str(), regex);
+}
+//! \endcond
+
+} // namespace gmx
diff --git a/src/gromacs/utility/gmxregex.h b/src/gromacs/utility/gmxregex.h
new file mode 100644 (file)
index 0000000..8de77c1
--- /dev/null
@@ -0,0 +1,131 @@
+/*
+ *
+ *                This source code is part of
+ *
+ *                 G   R   O   M   A   C   S
+ *
+ *          GROningen MAchine for Chemical Simulations
+ *
+ * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
+ * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
+ * Copyright (c) 2001-2009, The GROMACS development team,
+ * check out http://www.gromacs.org for more information.
+ *
+ * This program is free software; you can redistribute it and/or
+ * modify it under the terms of the GNU General Public License
+ * as published by the Free Software Foundation; either version 2
+ * of the License, or (at your option) any later version.
+ *
+ * If you want to redistribute modifications, 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 www.gromacs.org.
+ *
+ * To help us fund GROMACS development, we humbly ask that you cite
+ * the papers on the package - you can find them in the top README file.
+ *
+ * For more info, check our website at http://www.gromacs.org
+ */
+/*! \libinternal \file
+ * \brief
+ * Declares simple wrapper for regular expression functionality.
+ *
+ * \author Teemu Murtola <teemu.murtola@cbr.su.se>
+ * \inlibraryapi
+ * \ingroup module_utility
+ */
+#ifndef GMX_UTILITY_GMXREGEX_H
+#define GMX_UTILITY_GMXREGEX_H
+
+#include <string>
+
+#include "common.h"
+
+namespace gmx
+{
+
+class Regex;
+
+/*! \cond libapi */
+/*! \libinternal \brief
+ * Matches a string with a regular expression.
+ *
+ * \param[in] str   String to match.
+ * \param[in] regex Regular expression to match.
+ * \returns   true if \p regex matches the whole \p str.
+ *
+ * Does not throw currently, but this is subject to change if/when better error
+ * handling is implemented (currently, it returns false if the matching fails,
+ * e.g., because of out-of-memory).
+ *
+ * \ingroup module_utility
+ */
+bool regexMatch(const char *str, const Regex &regex);
+//! \copydoc regexMatch(const char *, const Regex &)
+bool regexMatch(const std::string &str, const Regex &regex);
+//! \endcond
+
+/*! \libinternal \brief
+ * Represents a regular expression.
+ *
+ * This class provides a simple interface for regular expression construction.
+ * regexMatch() is used to match the regular expression against a string.
+ * POSIX extended regular expression syntax is used.
+ *
+ * Currently, isSupported() will return false if POSIX regular expression
+ * header is not available (i.e., on Windows).  In this case, calling other
+ * constructors than the default constructor throws an exception.
+ *
+ * \see regexMatch()
+ *
+ * \inlibraryapi
+ * \ingroup module_utility
+ */
+class Regex
+{
+    public:
+        /*! \brief
+         * Returns true if regular expression support has been compiled in.
+         *
+         * Does not throw.
+         */
+        static bool isSupported();
+
+        /*! \brief
+         * Constructs a regular expression that matches nothing.
+         *
+         * Does not throw.
+         */
+        Regex();
+        /*! \brief
+         * Constructs a regular expression from a string.
+         *
+         * \param[in] value  String to compile into a regular expression.
+         * \throws    std::bad_alloc if out of memory.
+         * \throws    InvalidInputError if \p value is not a valid regular
+         *      expression.
+         *
+         * \todo
+         * Consider whether some other exception type would be better.
+         */
+        explicit Regex(const char *value);
+        //! \copydoc Regex(const char *)
+        explicit Regex(const std::string &value);
+        //! Frees memory allocated for the regular expression.
+        ~Regex();
+
+    private:
+        class Impl;
+
+        PrivateImplPointer<Impl> impl_;
+
+        friend bool regexMatch(const char *str, const Regex &regex);
+        friend bool regexMatch(const std::string &str, const Regex &regex);
+};
+
+} // namespace gmx
+
+#endif
+