3 * This source code is part of
7 * GROningen MAchine for Chemical Simulations
9 * Written by David van der Spoel, Erik Lindahl, Berk Hess, and others.
10 * Copyright (c) 1991-2000, University of Groningen, The Netherlands.
11 * Copyright (c) 2001-2009, The GROMACS development team,
12 * check out http://www.gromacs.org for more information.
14 * This program is free software; you can redistribute it and/or
15 * modify it under the terms of the GNU General Public License
16 * as published by the Free Software Foundation; either version 2
17 * of the License, or (at your option) any later version.
19 * If you want to redistribute modifications, please consider that
20 * scientific software is very special. Version control is crucial -
21 * bugs must be traceable. We will be happy to consider code for
22 * inclusion in the official distribution, but derived work must not
23 * be called official GROMACS. Details are found in the README & COPYING
24 * files - if they are missing, get the official version at www.gromacs.org.
26 * To help us fund GROMACS development, we humbly ask that you cite
27 * the papers on the package - you can find them in the top README file.
29 * For more info, check our website at http://www.gromacs.org
33 * Implements internal selection methods for numeric and string keyword
36 * \author Teemu Murtola <teemu.murtola@cbr.su.se>
37 * \ingroup module_selection
44 #include <boost/shared_ptr.hpp>
46 #include "gromacs/legacyheaders/macros.h"
47 #include "gromacs/legacyheaders/smalloc.h"
48 #include "gromacs/legacyheaders/string2.h"
50 #include "gromacs/selection/selmethod.h"
51 #include "gromacs/utility/exceptions.h"
52 #include "gromacs/utility/gmxregex.h"
53 #include "gromacs/utility/messagestringcollector.h"
54 #include "gromacs/utility/stringutil.h"
57 #include "parsetree.h"
61 /** Allocates data for integer keyword evaluation. */
63 init_data_kwint(int npar, gmx_ana_selparam_t *param);
64 /** Allocates data for real keyword evaluation. */
66 init_data_kwreal(int npar, gmx_ana_selparam_t *param);
67 /** Allocates data for string keyword evaluation. */
69 init_data_kwstr(int npar, gmx_ana_selparam_t *param);
70 /** Initializes data for integer keyword evaluation. */
72 init_kwint(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data);
73 /** Initializes data for real keyword evaluation. */
75 init_kwreal(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data);
76 /** Initializes data for string keyword evaluation. */
78 init_kwstr(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data);
79 /** Frees the memory allocated for string keyword evaluation. */
81 free_data_kwstr(void *data);
82 /** Evaluates integer selection keywords. */
84 evaluate_keyword_int(t_topology *top, t_trxframe *fr, t_pbc *pbc,
85 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data);
86 /** Evaluates real selection keywords. */
88 evaluate_keyword_real(t_topology *top, t_trxframe *fr, t_pbc *pbc,
89 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data);
90 /** Evaluates string selection keywords. */
92 evaluate_keyword_str(t_topology *top, t_trxframe *fr, t_pbc *pbc,
93 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data);
96 * Data structure for integer keyword expression evaluation.
98 typedef struct t_methoddata_kwint
100 /** Array of values for the keyword. */
102 /** Number of ranges in the \p r array. */
105 * Array of sorted integer ranges to match against.
107 * Each range is made of two integers, giving the endpoints (inclusive).
108 * This field stores the pointer to the ranges allocated by the
109 * parameter parser; see \ref SPAR_RANGES for more information.
112 } t_methoddata_kwint;
115 * Data structure for real keyword expression evaluation.
117 typedef struct t_methoddata_kwreal
119 /** Array of values for the keyword. */
121 /** Number of ranges in the \p r array. */
124 * Array of sorted ranges to match against.
126 * Each range is made of two values, giving the endpoints (inclusive).
127 * This field stores the pointer to the ranges allocated by the
128 * parameter parser; see \ref SPAR_RANGES for more information.
131 } t_methoddata_kwreal;
137 * Single item in the list of strings/regular expressions to match.
139 * \ingroup module_selection
141 class StringKeywordMatchItem
145 * Constructs a matcher from a string.
147 * \param[in] matchType String matching type.
148 * \param[in] str String to use for matching.
150 StringKeywordMatchItem(gmx::SelectionStringMatchType matchType,
154 bool bRegExp = (matchType == gmx::eStringMatchType_RegularExpression);
155 if (matchType == gmx::eStringMatchType_Auto)
157 for (size_t j = 0; j < std::strlen(str); ++j)
159 if (std::ispunct(str[j]) && str[j] != '?' && str[j] != '*')
168 if (!gmx::Regex::isSupported())
170 GMX_THROW(gmx::InvalidInputError(gmx::formatString(
171 "No regular expression support, "
172 "cannot match \"%s\"", str)));
174 regex_.reset(new gmx::Regex(str));
179 * Checks whether this item matches a string.
181 * \param[in] matchType String matching type.
182 * \param[in] value String to match.
183 * \returns true if this item matches \p value.
185 bool match(gmx::SelectionStringMatchType matchType,
186 const char *value) const
188 if (matchType == gmx::eStringMatchType_Exact)
190 return str_ == value;
194 return gmx::regexMatch(value, *regex_);
198 return gmx_wcmatch(str_.c_str(), value) == 0;
203 //! The raw string passed for the matcher.
205 //! Regular expression compiled from \p str_, if applicable.
206 boost::shared_ptr<gmx::Regex> regex_;
210 * Data structure for string keyword expression evaluation.
212 struct t_methoddata_kwstr
214 /** Matching type for the strings. */
215 gmx::SelectionStringMatchType matchType;
216 /** Array of values for the keyword. */
218 /** Array of strings/regular expressions to match against.*/
219 std::vector<StringKeywordMatchItem> matches;
224 /** Parameters for integer keyword evaluation. */
225 static gmx_ana_selparam_t smparams_keyword_int[] = {
226 {NULL, {INT_VALUE, -1, {NULL}}, NULL, SPAR_ATOMVAL},
227 {NULL, {INT_VALUE, -1, {NULL}}, NULL, SPAR_RANGES | SPAR_VARNUM},
230 /** Parameters for real keyword evaluation. */
231 static gmx_ana_selparam_t smparams_keyword_real[] = {
232 {NULL, {REAL_VALUE, -1, {NULL}}, NULL, SPAR_ATOMVAL | SPAR_DYNAMIC},
233 {NULL, {REAL_VALUE, -1, {NULL}}, NULL, SPAR_RANGES | SPAR_VARNUM},
236 /** Parameters for string keyword evaluation. */
237 static gmx_ana_selparam_t smparams_keyword_str[] = {
238 {NULL, {STR_VALUE, -1, {NULL}}, NULL, SPAR_ATOMVAL},
239 {NULL, {STR_VALUE, -1, {NULL}}, NULL, SPAR_VARNUM},
242 /** \internal Selection method data for integer keyword evaluation. */
243 gmx_ana_selmethod_t sm_keyword_int = {
244 "kw_int", GROUP_VALUE, SMETH_SINGLEVAL,
245 asize(smparams_keyword_int), smparams_keyword_int,
252 &evaluate_keyword_int,
257 /** \internal Selection method data for real keyword evaluation. */
258 gmx_ana_selmethod_t sm_keyword_real = {
259 "kw_real", GROUP_VALUE, SMETH_SINGLEVAL,
260 asize(smparams_keyword_real), smparams_keyword_real,
267 &evaluate_keyword_real,
272 /** \internal Selection method data for string keyword evaluation. */
273 gmx_ana_selmethod_t sm_keyword_str = {
274 "kw_str", GROUP_VALUE, SMETH_SINGLEVAL,
275 asize(smparams_keyword_str), smparams_keyword_str,
282 &evaluate_keyword_str,
287 /** Initializes keyword evaluation for an arbitrary group. */
289 init_kweval(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data);
290 /** Initializes output for keyword evaluation in an arbitrary group. */
292 init_output_kweval(t_topology *top, gmx_ana_selvalue_t *out, void *data);
293 /** Frees the data allocated for keyword evaluation in an arbitrary group. */
295 free_data_kweval(void *data);
296 /** Initializes frame evaluation for keyword evaluation in an arbitrary group. */
298 init_frame_kweval(t_topology *top, t_trxframe *fr, t_pbc *pbc, void *data);
299 /** Evaluates keywords in an arbitrary group. */
301 evaluate_kweval(t_topology *top, t_trxframe *fr, t_pbc *pbc,
302 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data);
305 * Data structure for keyword evaluation in arbitrary groups.
309 /** Wrapped keyword method for evaluating the values. */
310 gmx_ana_selmethod_t *kwmethod;
311 /** Method data for \p kwmethod. */
313 /** Group in which \p kwmethod should be evaluated. */
315 } t_methoddata_kweval;
317 /** Parameters for keyword evaluation in an arbitrary group. */
318 static gmx_ana_selparam_t smparams_kweval[] = {
319 {NULL, {GROUP_VALUE, 1, {NULL}}, NULL, SPAR_DYNAMIC},
323 /********************************************************************
324 * INTEGER KEYWORD EVALUATION
325 ********************************************************************/
328 * \param[in] npar Not used.
329 * \param param Not used.
330 * \returns Pointer to the allocated data (\ref t_methoddata_kwint).
332 * Allocates memory for a \ref t_methoddata_kwint structure.
335 init_data_kwint(int npar, gmx_ana_selparam_t *param)
337 t_methoddata_kwint *data;
344 * \param[in] top Not used.
345 * \param[in] npar Not used (should be 2).
346 * \param[in] param Method parameters (should point to \ref smparams_keyword_int).
347 * \param[in] data Should point to \ref t_methoddata_kwint.
350 init_kwint(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data)
352 t_methoddata_kwint *d = (t_methoddata_kwint *)data;
354 d->v = param[0].val.u.i;
355 d->n = param[1].val.nr;
356 d->r = param[1].val.u.i;
360 * See sel_updatefunc() for description of the parameters.
361 * \p data should point to a \c t_methoddata_kwint.
363 * Does a binary search to find which atoms match the ranges in the
364 * \c t_methoddata_kwint structure for this selection.
365 * Matching atoms are stored in \p out->u.g.
368 evaluate_keyword_int(t_topology *top, t_trxframe *fr, t_pbc *pbc,
369 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
371 t_methoddata_kwint *d = (t_methoddata_kwint *)data;
372 int n, i, j, jmin, jmax;
377 for (i = 0; i < g->isize; ++i)
380 if (d->r[0] > val || d->r[2*n-1] < val)
386 while (jmax - jmin > 1)
388 j = jmin + (jmax - jmin) / 2;
396 if (val <= d->r[2*j+1])
403 if (val <= d->r[2*jmin+1])
405 out->u.g->index[out->u.g->isize++] = g->index[i];
411 /********************************************************************
412 * REAL KEYWORD EVALUATION
413 ********************************************************************/
416 * \param[in] npar Not used.
417 * \param param Not used.
418 * \returns Pointer to the allocated data (\ref t_methoddata_kwreal).
420 * Allocates memory for a \ref t_methoddata_kwreal structure.
423 init_data_kwreal(int npar, gmx_ana_selparam_t *param)
425 t_methoddata_kwreal *data;
432 * \param[in] top Not used.
433 * \param[in] npar Not used (should be 2).
434 * \param[in] param Method parameters (should point to \ref smparams_keyword_real).
435 * \param[in] data Should point to \ref t_methoddata_kwreal.
436 * \returns 0 (the initialization always succeeds).
439 init_kwreal(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data)
441 t_methoddata_kwreal *d = (t_methoddata_kwreal *)data;
443 d->v = param[0].val.u.r;
444 d->n = param[1].val.nr;
445 d->r = param[1].val.u.r;
449 * See sel_updatefunc() for description of the parameters.
450 * \p data should point to a \c t_methoddata_kwreal.
452 * Does a binary search to find which atoms match the ranges in the
453 * \c t_methoddata_kwreal structure for this selection.
454 * Matching atoms are stored in \p out->u.g.
457 evaluate_keyword_real(t_topology *top, t_trxframe *fr, t_pbc *pbc,
458 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
460 t_methoddata_kwreal *d = (t_methoddata_kwreal *)data;
461 int n, i, j, jmin, jmax;
466 for (i = 0; i < g->isize; ++i)
469 if (d->r[0] > val || d->r[2*n-1] < val)
475 while (jmax - jmin > 1)
477 j = jmin + (jmax - jmin) / 2;
485 if (val <= d->r[2*j+1])
492 if (val <= d->r[2*jmin+1])
494 out->u.g->index[out->u.g->isize++] = g->index[i];
500 /********************************************************************
501 * STRING KEYWORD EVALUATION
502 ********************************************************************/
505 * \param[in] npar Not used.
506 * \param param Not used.
507 * \returns Pointer to the allocated data (t_methoddata_kwstr).
509 * Allocates memory for a t_methoddata_kwstr structure.
512 init_data_kwstr(int npar, gmx_ana_selparam_t *param)
514 t_methoddata_kwstr *data = new t_methoddata_kwstr();
515 data->matchType = gmx::eStringMatchType_Auto;
520 * \param[in,out] sel Selection element to initialize.
521 * \param[in] matchType Method to use to match string values.
523 * Sets the string matching method for string keyword matching.
526 _gmx_selelem_set_kwstr_match_type(const gmx::SelectionTreeElementPointer &sel,
527 gmx::SelectionStringMatchType matchType)
529 t_methoddata_kwstr *d = (t_methoddata_kwstr *)sel->u.expr.mdata;
531 if (sel->type != SEL_EXPRESSION || !sel->u.expr.method
532 || sel->u.expr.method->name != sm_keyword_str.name)
536 d->matchType = matchType;
540 * \param[in] top Not used.
541 * \param[in] npar Not used (should be 2).
542 * \param[in] param Method parameters (should point to \ref smparams_keyword_str).
543 * \param[in] data Should point to t_methoddata_kwstr.
546 init_kwstr(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data)
548 t_methoddata_kwstr *d = (t_methoddata_kwstr *)data;
550 d->v = param[0].val.u.s;
551 /* Return if this is not the first time */
552 if (!d->matches.empty())
556 int n = param[1].val.nr;
557 d->matches.reserve(n);
558 for (int i = 0; i < n; ++i)
560 const char *s = param[1].val.u.s[i];
561 d->matches.push_back(StringKeywordMatchItem(d->matchType, s));
566 * \param data Data to free (should point to a t_methoddata_kwstr).
569 free_data_kwstr(void *data)
571 t_methoddata_kwstr *d = (t_methoddata_kwstr *)data;
576 * See sel_updatefunc() for description of the parameters.
577 * \p data should point to a \c t_methoddata_kwstr.
579 * Does a linear search to find which atoms match the strings in the
580 * \c t_methoddata_kwstr structure for this selection.
581 * Wildcards are allowed in the strings.
582 * Matching atoms are stored in \p out->u.g.
585 evaluate_keyword_str(t_topology *top, t_trxframe *fr, t_pbc *pbc,
586 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
588 t_methoddata_kwstr *d = (t_methoddata_kwstr *)data;
591 for (int i = 0; i < g->isize; ++i)
593 for (size_t j = 0; j < d->matches.size(); ++j)
595 if (d->matches[j].match(d->matchType, d->v[i]))
597 out->u.g->index[out->u.g->isize++] = g->index[i];
605 /********************************************************************
606 * KEYWORD EVALUATION FOR ARBITRARY GROUPS
607 ********************************************************************/
610 * \param[in] top Not used.
611 * \param[in] npar Not used.
612 * \param[in] param Not used.
613 * \param[in] data Should point to \ref t_methoddata_kweval.
614 * \returns 0 on success, a non-zero error code on return.
616 * Calls the initialization method of the wrapped keyword.
619 init_kweval(t_topology *top, int npar, gmx_ana_selparam_t *param, void *data)
621 t_methoddata_kweval *d = (t_methoddata_kweval *)data;
623 d->kwmethod->init(top, 0, NULL, d->kwmdata);
627 * \param[in] top Not used.
628 * \param[in,out] out Pointer to output data structure.
629 * \param[in,out] data Should point to \c t_methoddata_kweval.
630 * \returns 0 for success.
633 init_output_kweval(t_topology *top, gmx_ana_selvalue_t *out, void *data)
635 t_methoddata_kweval *d = (t_methoddata_kweval *)data;
637 out->nr = d->g.isize;
638 _gmx_selvalue_reserve(out, out->nr);
642 * \param data Data to free (should point to a \c t_methoddata_kweval).
644 * Frees the memory allocated for all the members of \c t_methoddata_kweval.
647 free_data_kweval(void *data)
649 t_methoddata_kweval *d = (t_methoddata_kweval *)data;
651 _gmx_selelem_free_method(d->kwmethod, d->kwmdata);
656 * \param[in] top Topology.
657 * \param[in] fr Current frame.
658 * \param[in] pbc PBC structure.
659 * \param data Should point to a \ref t_methoddata_kweval.
660 * \returns 0 on success, a non-zero error code on error.
662 * Creates a lookup structure that enables fast queries of whether a point
663 * is within the solid angle or not.
666 init_frame_kweval(t_topology *top, t_trxframe *fr, t_pbc *pbc, void *data)
668 t_methoddata_kweval *d = (t_methoddata_kweval *)data;
670 d->kwmethod->init_frame(top, fr, pbc, d->kwmdata);
674 * See sel_updatefunc() for description of the parameters.
675 * \p data should point to a \c t_methoddata_kweval.
677 * Calls the evaluation function of the wrapped keyword with the given
678 * parameters, with the exception of using \c t_methoddata_kweval::g for the
682 evaluate_kweval(t_topology *top, t_trxframe *fr, t_pbc *pbc,
683 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
685 t_methoddata_kweval *d = (t_methoddata_kweval *)data;
687 d->kwmethod->update(top, fr, pbc, &d->g, out, d->kwmdata);
691 * \param[in] method Keyword selection method to evaluate.
692 * \param[in] params Parameter that gives the group to evaluate \p method in.
693 * \param[in] scanner Scanner data structure.
694 * \returns Pointer to the created selection element (NULL on error).
696 * Creates a \ref SEL_EXPRESSION selection element that evaluates the keyword
697 * method given by \p method in the group given by \p param.
699 * The name of \p param should be empty.
701 gmx::SelectionTreeElementPointer
702 _gmx_sel_init_keyword_evaluator(gmx_ana_selmethod_t *method,
703 const gmx::SelectionParserParameterList ¶ms,
706 gmx::MessageStringCollector *errors = _gmx_sel_lexer_error_reporter(scanner);
708 sprintf(buf, "In evaluation of '%s'", method->name);
709 gmx::MessageStringContext context(errors, buf);
711 if ((method->flags & (SMETH_SINGLEVAL | SMETH_VARNUMVAL))
712 || method->outinit || method->pupdate)
714 GMX_THROW(gmx::InternalError(
715 "Unsupported keyword method for arbitrary group evaluation"));
718 gmx::SelectionTreeElementPointer sel(
719 new gmx::SelectionTreeElement(SEL_EXPRESSION));
720 _gmx_selelem_set_method(sel, method, scanner);
722 t_methoddata_kweval *data;
724 data->kwmethod = sel->u.expr.method;
725 data->kwmdata = sel->u.expr.mdata;
726 gmx_ana_index_clear(&data->g);
728 snew(sel->u.expr.method, 1);
729 memcpy(sel->u.expr.method, data->kwmethod, sizeof(gmx_ana_selmethod_t));
730 sel->u.expr.method->flags |= SMETH_VARNUMVAL;
731 sel->u.expr.method->init_data = NULL;
732 sel->u.expr.method->set_poscoll = NULL;
733 sel->u.expr.method->init = method->init ? &init_kweval : NULL;
734 sel->u.expr.method->outinit = &init_output_kweval;
735 sel->u.expr.method->free = &free_data_kweval;
736 sel->u.expr.method->init_frame = method->init_frame ? &init_frame_kweval : NULL;
737 sel->u.expr.method->update = &evaluate_kweval;
738 sel->u.expr.method->pupdate = NULL;
739 sel->u.expr.method->nparams = asize(smparams_kweval);
740 sel->u.expr.method->param = smparams_kweval;
741 _gmx_selelem_init_method_params(sel, scanner);
742 sel->u.expr.mdata = data;
744 sel->u.expr.method->param[0].val.u.g = &data->g;
746 if (!_gmx_sel_parse_params(params, sel->u.expr.method->nparams,
747 sel->u.expr.method->param, sel, scanner))
749 return gmx::SelectionTreeElementPointer();