2 * This file is part of the GROMACS molecular simulation package.
4 * Copyright (c) 2009,2010,2011,2012,2013,2014,2015,2016, by the GROMACS development team, led by
5 * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
6 * and including many others, as listed in the AUTHORS file in the
7 * top-level source directory and at http://www.gromacs.org.
9 * GROMACS is free software; you can redistribute it and/or
10 * modify it under the terms of the GNU Lesser General Public License
11 * as published by the Free Software Foundation; either version 2.1
12 * of the License, or (at your option) any later version.
14 * GROMACS is distributed in the hope that it will be useful,
15 * but WITHOUT ANY WARRANTY; without even the implied warranty of
16 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
17 * Lesser General Public License for more details.
19 * You should have received a copy of the GNU Lesser General Public
20 * License along with GROMACS; if not, see
21 * http://www.gnu.org/licenses, or write to the Free Software Foundation,
22 * Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA.
24 * If you want to redistribute modifications to GROMACS, please
25 * consider that scientific software is very special. Version
26 * control is crucial - bugs must be traceable. We will be happy to
27 * consider code for inclusion in the official distribution, but
28 * derived work must not be called official GROMACS. Details are found
29 * in the README & COPYING files - if they are missing, get the
30 * official version at http://www.gromacs.org.
32 * To help us fund GROMACS development, we humbly ask that you cite
33 * the research papers on the package. Check out http://www.gromacs.org.
37 * Implements internal selection methods for numeric and string keyword
40 * \author Teemu Murtola <teemu.murtola@gmail.com>
41 * \ingroup module_selection
51 #include "gromacs/selection/position.h"
52 #include "gromacs/utility/arraysize.h"
53 #include "gromacs/utility/cstringutil.h"
54 #include "gromacs/utility/exceptions.h"
55 #include "gromacs/utility/gmxassert.h"
56 #include "gromacs/utility/gmxregex.h"
57 #include "gromacs/utility/smalloc.h"
58 #include "gromacs/utility/stringutil.h"
61 #include "parsetree.h"
64 #include "selmethod.h"
67 * Allocates data for integer keyword evaluation.
69 * \param[in] npar Not used.
70 * \param param Not used.
71 * \returns Pointer to the allocated data (\ref t_methoddata_kwint).
73 * Allocates memory for a \ref t_methoddata_kwint structure.
76 init_data_kwint(int npar, gmx_ana_selparam_t * param);
78 * Allocates data for real keyword evaluation.
80 * \param[in] npar Not used.
81 * \param param Not used.
82 * \returns Pointer to the allocated data (\ref t_methoddata_kwreal).
84 * Allocates memory for a \ref t_methoddata_kwreal structure.
87 init_data_kwreal(int npar, gmx_ana_selparam_t * param);
89 * Allocates data for string keyword evaluation.
91 * \param[in] npar Not used.
92 * \param param Not used.
93 * \returns Pointer to the allocated data (t_methoddata_kwstr).
95 * Allocates memory for a t_methoddata_kwstr structure.
98 init_data_kwstr(int npar, gmx_ana_selparam_t * param);
99 /** /brief Initializes data for integer keyword evaluation.
101 * \param[in] top Not used.
102 * \param[in] npar Not used (should be 2).
103 * \param[in] param Method parameters (should point to \ref smparams_keyword_int).
104 * \param[in] data Should point to \ref t_methoddata_kwint.
107 init_kwint(const gmx_mtop_t *top, int npar, gmx_ana_selparam_t *param, void *data);
109 * Initializes data for real keyword evaluation.
111 * \param[in] top Not used.
112 * \param[in] npar Not used (should be 2).
113 * \param[in] param Method parameters (should point to \ref smparams_keyword_real).
114 * \param[in] data Should point to \ref t_methoddata_kwreal.
115 * \returns 0 (the initialization always succeeds).
118 init_kwreal(const gmx_mtop_t *top, int npar, gmx_ana_selparam_t *param, void *data);
120 * Initializes data for string keyword evaluation.
122 * \param[in] top Not used.
123 * \param[in] npar Not used (should be 2).
124 * \param[in] param Method parameters (should point to \ref smparams_keyword_str).
125 * \param[in] data Should point to t_methoddata_kwstr.
128 init_kwstr(const gmx_mtop_t *top, int npar, gmx_ana_selparam_t *param, void *data);
129 /** Frees the memory allocated for string keyword evaluation. */
131 free_data_kwstr(void *data);
132 /** Evaluates integer selection keywords. */
134 evaluate_keyword_int(const gmx::SelMethodEvalContext &context,
135 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data);
136 /** Evaluates real selection keywords. */
138 evaluate_keyword_real(const gmx::SelMethodEvalContext &context,
139 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data);
140 /** Evaluates string selection keywords. */
142 evaluate_keyword_str(const gmx::SelMethodEvalContext &context,
143 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data);
146 * Data structure for integer keyword expression evaluation.
148 typedef struct t_methoddata_kwint
150 /** Array of values for the keyword. */
152 /** Number of ranges in the \p r array. */
155 * Array of sorted integer ranges to match against.
157 * Each range is made of two integers, giving the endpoints (inclusive).
158 * This field stores the pointer to the ranges allocated by the
159 * parameter parser; see \ref SPAR_RANGES for more information.
162 } t_methoddata_kwint;
165 * Data structure for real keyword expression evaluation.
167 typedef struct t_methoddata_kwreal
169 /** Array of values for the keyword. */
171 /** Number of ranges in the \p r array. */
174 * Array of sorted ranges to match against.
176 * Each range is made of two values, giving the endpoints (inclusive).
177 * This field stores the pointer to the ranges allocated by the
178 * parameter parser; see \ref SPAR_RANGES for more information.
181 } t_methoddata_kwreal;
187 * Single item in the list of strings/regular expressions to match.
189 * \ingroup module_selection
191 class StringKeywordMatchItem
195 * Constructs a matcher from a string.
197 * \param[in] matchType String matching type.
198 * \param[in] str String to use for matching.
200 StringKeywordMatchItem(gmx::SelectionStringMatchType matchType,
204 bool bRegExp = (matchType == gmx::eStringMatchType_RegularExpression);
205 if (matchType == gmx::eStringMatchType_Auto)
207 for (size_t j = 0; j < std::strlen(str); ++j)
209 if (std::ispunct(str[j]) && str[j] != '?' && str[j] != '*')
218 if (!gmx::Regex::isSupported())
221 = gmx::formatString("No regular expression support, "
222 "cannot match \"%s\"", str);
223 GMX_THROW(gmx::InvalidInputError(message));
225 regex_.reset(new gmx::Regex(str));
230 * Checks whether this item matches a string.
232 * \param[in] matchType String matching type.
233 * \param[in] value String to match.
234 * \returns true if this item matches \p value.
236 bool match(gmx::SelectionStringMatchType matchType,
237 const char *value) const
239 if (matchType == gmx::eStringMatchType_Exact)
241 return str_ == value;
245 return gmx::regexMatch(value, *regex_);
249 return gmx_wcmatch(str_.c_str(), value) == 0;
254 //! The raw string passed for the matcher.
256 //! Regular expression compiled from \p str_, if applicable.
257 std::shared_ptr<gmx::Regex> regex_;
261 * Data structure for string keyword expression evaluation.
263 struct t_methoddata_kwstr
265 /** Matching type for the strings. */
266 gmx::SelectionStringMatchType matchType;
267 /** Array of values for the keyword. */
269 /** Array of strings/regular expressions to match against.*/
270 std::vector<StringKeywordMatchItem> matches;
275 /** Parameters for integer keyword evaluation. */
276 static gmx_ana_selparam_t smparams_keyword_int[] = {
277 {NULL, {INT_VALUE, -1, {NULL}}, NULL, SPAR_ATOMVAL},
278 {NULL, {INT_VALUE, -1, {NULL}}, NULL, SPAR_RANGES | SPAR_VARNUM},
281 /** Parameters for real keyword evaluation. */
282 static gmx_ana_selparam_t smparams_keyword_real[] = {
283 {NULL, {REAL_VALUE, -1, {NULL}}, NULL, SPAR_ATOMVAL | SPAR_DYNAMIC},
284 {NULL, {REAL_VALUE, -1, {NULL}}, NULL, SPAR_RANGES | SPAR_VARNUM},
287 /** Parameters for string keyword evaluation. */
288 static gmx_ana_selparam_t smparams_keyword_str[] = {
289 {NULL, {STR_VALUE, -1, {NULL}}, NULL, SPAR_ATOMVAL},
290 {NULL, {STR_VALUE, -1, {NULL}}, NULL, SPAR_VARNUM},
293 /** Selection method data for integer keyword evaluation. */
294 gmx_ana_selmethod_t sm_keyword_int = {
295 "kw_int", GROUP_VALUE, SMETH_SINGLEVAL,
296 asize(smparams_keyword_int), smparams_keyword_int,
303 &evaluate_keyword_int,
305 {NULL, NULL, 0, NULL},
308 /** Selection method data for real keyword evaluation. */
309 gmx_ana_selmethod_t sm_keyword_real = {
310 "kw_real", GROUP_VALUE, SMETH_SINGLEVAL,
311 asize(smparams_keyword_real), smparams_keyword_real,
318 &evaluate_keyword_real,
320 {NULL, NULL, 0, NULL},
323 /** Selection method data for string keyword evaluation. */
324 gmx_ana_selmethod_t sm_keyword_str = {
325 "kw_str", GROUP_VALUE, SMETH_SINGLEVAL,
326 asize(smparams_keyword_str), smparams_keyword_str,
333 &evaluate_keyword_str,
335 {NULL, NULL, 0, NULL},
339 * Initializes keyword evaluation for an arbitrary group.
341 * \param[in] top Not used.
342 * \param[in] npar Not used.
343 * \param[in] param Not used.
344 * \param[in] data Should point to \ref t_methoddata_kweval.
345 * \returns 0 on success, a non-zero error code on return.
347 * Calls the initialization method of the wrapped keyword.
350 init_kweval(const gmx_mtop_t *top, int npar, gmx_ana_selparam_t * param, void *data);
352 * Initializes output for keyword evaluation in an arbitrary group.
354 * \param[in] top Not used.
355 * \param[in,out] out Pointer to output data structure.
356 * \param[in,out] data Should point to \c t_methoddata_kweval.
357 * \returns 0 for success.
360 init_output_kweval(const gmx_mtop_t *top, gmx_ana_selvalue_t *out, void *data);
361 /** Frees the data allocated for keyword evaluation in an arbitrary group. */
363 free_data_kweval(void *data);
364 /** Initializes frame evaluation for keyword evaluation in an arbitrary group. */
366 init_frame_kweval(const gmx::SelMethodEvalContext &context, void *data);
368 * Evaluates keywords in an arbitrary group.
370 * See sel_updatefunc() for description of the parameters.
371 * \p data should point to a \c t_methoddata_kweval.
373 * Calls the evaluation function of the wrapped keyword with the given
374 * parameters, with the exception of using \c t_methoddata_kweval::g for the
378 evaluate_kweval(const gmx::SelMethodEvalContext &context, gmx_ana_index_t *g,
379 gmx_ana_selvalue_t *out, void *data);
381 * Evaluates keywords in an arbitrary set of positions.
383 * See sel_updatefunc() for description of the parameters.
384 * \p data should point to a \c t_methoddata_kweval.
386 * Calls the evaluation function of the wrapped keyword with the given
387 * parameters, with the exception of using \c t_methoddata_kweval::p for the
388 * evaluation positions.
391 evaluate_kweval_pos(const gmx::SelMethodEvalContext &context, gmx_ana_index_t *g,
392 gmx_ana_selvalue_t *out, void *data);
395 * Data structure for keyword evaluation in arbitrary groups.
397 struct t_methoddata_kweval
399 //! Initialize new keyword evaluator for the given keyword.
400 t_methoddata_kweval(gmx_ana_selmethod_t *method, void *data)
401 : kwmethod(method), kwmdata(data)
403 gmx_ana_index_clear(&g);
405 ~t_methoddata_kweval()
407 _gmx_selelem_free_method(kwmethod, kwmdata);
410 //! Wrapped keyword method for evaluating the values.
411 gmx_ana_selmethod_t *kwmethod;
412 //! Method data for \p kwmethod.
414 //! Group in which \p kwmethod should be evaluated.
416 //! Positions for which \p kwmethod should be evaluated.
420 /** Parameters for keyword evaluation in an arbitrary group. */
421 static gmx_ana_selparam_t smparams_kweval_group[] = {
422 {NULL, {GROUP_VALUE, 1, {NULL}}, NULL, SPAR_DYNAMIC},
424 /** Parameters for keyword evaluation from positions. */
425 static gmx_ana_selparam_t smparams_kweval_pos[] = {
426 {NULL, {POS_VALUE, 1, {NULL}}, NULL, SPAR_DYNAMIC},
430 /********************************************************************
431 * INTEGER KEYWORD EVALUATION
432 ********************************************************************/
435 init_data_kwint(int /* npar */, gmx_ana_selparam_t * /* param */)
437 t_methoddata_kwint *data;
444 init_kwint(const gmx_mtop_t * /* top */, int /* npar */, gmx_ana_selparam_t *param, void *data)
446 t_methoddata_kwint *d = (t_methoddata_kwint *)data;
448 d->v = param[0].val.u.i;
449 d->n = param[1].val.nr;
450 d->r = param[1].val.u.i;
454 * See sel_updatefunc() for description of the parameters.
455 * \p data should point to a \c t_methoddata_kwint.
457 * Does a binary search to find which atoms match the ranges in the
458 * \c t_methoddata_kwint structure for this selection.
459 * Matching atoms are stored in \p out->u.g.
462 evaluate_keyword_int(const gmx::SelMethodEvalContext & /*context*/,
463 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
465 t_methoddata_kwint *d = (t_methoddata_kwint *)data;
466 int n, i, j, jmin, jmax;
471 for (i = 0; i < g->isize; ++i)
474 if (d->r[0] > val || d->r[2*n-1] < val)
480 while (jmax - jmin > 1)
482 j = jmin + (jmax - jmin) / 2;
490 if (val <= d->r[2*j+1])
497 if (val <= d->r[2*jmin+1])
499 out->u.g->index[out->u.g->isize++] = g->index[i];
505 /********************************************************************
506 * REAL KEYWORD EVALUATION
507 ********************************************************************/
510 init_data_kwreal(int /* npar */, gmx_ana_selparam_t * /* param */)
512 t_methoddata_kwreal *data;
519 init_kwreal(const gmx_mtop_t * /* top */, int /* npar */, gmx_ana_selparam_t *param, void *data)
521 t_methoddata_kwreal *d = (t_methoddata_kwreal *)data;
523 d->v = param[0].val.u.r;
524 d->n = param[1].val.nr;
525 d->r = param[1].val.u.r;
529 * See sel_updatefunc() for description of the parameters.
530 * \p data should point to a \c t_methoddata_kwreal.
532 * Does a binary search to find which atoms match the ranges in the
533 * \c t_methoddata_kwreal structure for this selection.
534 * Matching atoms are stored in \p out->u.g.
537 evaluate_keyword_real(const gmx::SelMethodEvalContext & /*context*/,
538 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
540 t_methoddata_kwreal *d = (t_methoddata_kwreal *)data;
541 int n, i, j, jmin, jmax;
546 for (i = 0; i < g->isize; ++i)
549 if (d->r[0] > val || d->r[2*n-1] < val)
555 while (jmax - jmin > 1)
557 j = jmin + (jmax - jmin) / 2;
565 if (val <= d->r[2*j+1])
572 if (val <= d->r[2*jmin+1])
574 out->u.g->index[out->u.g->isize++] = g->index[i];
580 /********************************************************************
581 * STRING KEYWORD EVALUATION
582 ********************************************************************/
585 init_data_kwstr(int /* npar */, gmx_ana_selparam_t * /* param */)
587 t_methoddata_kwstr *data = new t_methoddata_kwstr();
588 data->matchType = gmx::eStringMatchType_Auto;
593 * \param[in,out] sel Selection element to initialize.
594 * \param[in] matchType Method to use to match string values.
596 * Sets the string matching method for string keyword matching.
599 _gmx_selelem_set_kwstr_match_type(const gmx::SelectionTreeElementPointer &sel,
600 gmx::SelectionStringMatchType matchType)
602 t_methoddata_kwstr *d = static_cast<t_methoddata_kwstr *>(sel->u.expr.mdata);
604 if (sel->type != SEL_EXPRESSION || !sel->u.expr.method
605 || sel->u.expr.method->name != sm_keyword_str.name)
609 d->matchType = matchType;
613 init_kwstr(const gmx_mtop_t * /* top */, int /* npar */, gmx_ana_selparam_t *param, void *data)
615 t_methoddata_kwstr *d = static_cast<t_methoddata_kwstr *>(data);
617 d->v = param[0].val.u.s;
618 /* Return if this is not the first time */
619 if (!d->matches.empty())
623 int n = param[1].val.nr;
624 d->matches.reserve(n);
625 for (int i = 0; i < n; ++i)
627 const char *s = param[1].val.u.s[i];
628 d->matches.emplace_back(d->matchType, s);
633 * \param data Data to free (should point to a t_methoddata_kwstr).
636 free_data_kwstr(void *data)
638 t_methoddata_kwstr *d = static_cast<t_methoddata_kwstr *>(data);
643 * See sel_updatefunc() for description of the parameters.
644 * \p data should point to a \c t_methoddata_kwstr.
646 * Does a linear search to find which atoms match the strings in the
647 * \c t_methoddata_kwstr structure for this selection.
648 * Wildcards are allowed in the strings.
649 * Matching atoms are stored in \p out->u.g.
652 evaluate_keyword_str(const gmx::SelMethodEvalContext & /*context*/,
653 gmx_ana_index_t *g, gmx_ana_selvalue_t *out, void *data)
655 t_methoddata_kwstr *d = static_cast<t_methoddata_kwstr *>(data);
658 for (int i = 0; i < g->isize; ++i)
660 for (size_t j = 0; j < d->matches.size(); ++j)
662 if (d->matches[j].match(d->matchType, d->v[i]))
664 out->u.g->index[out->u.g->isize++] = g->index[i];
672 /********************************************************************
673 * KEYWORD EVALUATION FOR ARBITRARY GROUPS
674 ********************************************************************/
677 init_kweval(const gmx_mtop_t *top, int /* npar */, gmx_ana_selparam_t * /* param */, void *data)
679 t_methoddata_kweval *d = static_cast<t_methoddata_kweval *>(data);
681 d->kwmethod->init(top, 0, NULL, d->kwmdata);
685 init_output_kweval(const gmx_mtop_t * /* top */, gmx_ana_selvalue_t *out, void *data)
687 t_methoddata_kweval *d = static_cast<t_methoddata_kweval *>(data);
689 out->nr = d->g.isize;
690 _gmx_selvalue_reserve(out, out->nr);
694 * \param data Data to free (should point to a \c t_methoddata_kweval).
696 * Frees the memory allocated for all the members of \c t_methoddata_kweval.
699 free_data_kweval(void *data)
701 t_methoddata_kweval *d = static_cast<t_methoddata_kweval *>(data);
707 * \param[in] context Evaluation context.
708 * \param data Should point to a \ref t_methoddata_kweval.
711 init_frame_kweval(const gmx::SelMethodEvalContext &context, void *data)
713 t_methoddata_kweval *d = static_cast<t_methoddata_kweval *>(data);
715 d->kwmethod->init_frame(context, d->kwmdata);
719 evaluate_kweval(const gmx::SelMethodEvalContext &context,
720 gmx_ana_index_t * /* g */, gmx_ana_selvalue_t *out, void *data)
722 t_methoddata_kweval *d = static_cast<t_methoddata_kweval *>(data);
724 d->kwmethod->update(context, &d->g, out, d->kwmdata);
728 evaluate_kweval_pos(const gmx::SelMethodEvalContext &context,
729 gmx_ana_index_t * /* g */, gmx_ana_selvalue_t *out, void *data)
731 t_methoddata_kweval *d = static_cast<t_methoddata_kweval *>(data);
733 d->kwmethod->pupdate(context, &d->p, out, d->kwmdata);
737 * Initializes a selection element for evaluating a keyword in a given group.
739 * \param[in] method Keyword selection method to evaluate.
740 * \param[in] params Parameters to pass to initialization (the child group).
741 * \param[in] scanner Scanner data structure.
742 * \returns Pointer to the created selection element.
744 * Implements _gmx_sel_init_keyword_evaluator() for \ref GROUP_VALUE input
747 static gmx::SelectionTreeElementPointer
748 init_evaluator_group(gmx_ana_selmethod_t *method,
749 const gmx::SelectionParserParameterList ¶ms,
752 if ((method->flags & (SMETH_SINGLEVAL | SMETH_VARNUMVAL))
753 || method->outinit || method->pupdate)
756 = gmx::formatString("Keyword '%s' cannot be evaluated in this context",
758 GMX_THROW(gmx::InvalidInputError(message));
761 // TODO: For same ... as ..., some other location could be more intuitive.
762 gmx::SelectionTreeElementPointer sel(
763 new gmx::SelectionTreeElement(
764 SEL_EXPRESSION, _gmx_sel_lexer_get_current_location(scanner)));
765 _gmx_selelem_set_method(sel, method, scanner);
767 t_methoddata_kweval *data
768 = new t_methoddata_kweval(sel->u.expr.method, sel->u.expr.mdata);
770 snew(sel->u.expr.method, 1);
771 sel->u.expr.method->name = data->kwmethod->name;
772 sel->u.expr.method->type = data->kwmethod->type;
773 sel->u.expr.method->flags = data->kwmethod->flags | SMETH_VARNUMVAL;
774 sel->u.expr.method->init_data = NULL;
775 sel->u.expr.method->set_poscoll = NULL;
776 sel->u.expr.method->init = method->init ? &init_kweval : NULL;
777 sel->u.expr.method->outinit = &init_output_kweval;
778 sel->u.expr.method->free = &free_data_kweval;
779 sel->u.expr.method->init_frame = method->init_frame ? &init_frame_kweval : NULL;
780 sel->u.expr.method->update = &evaluate_kweval;
781 sel->u.expr.method->pupdate = NULL;
782 sel->u.expr.method->nparams = asize(smparams_kweval_group);
783 sel->u.expr.method->param = smparams_kweval_group;
784 _gmx_selelem_init_method_params(sel, scanner);
785 sel->u.expr.mdata = data;
787 sel->u.expr.method->param[0].val.u.g = &data->g;
789 _gmx_sel_parse_params(params, sel->u.expr.method->nparams,
790 sel->u.expr.method->param, sel, scanner);
795 * Initializes a selection element for evaluating a keyword in given positions.
797 * \param[in] method Keyword selection method to evaluate.
798 * \param[in] params Parameters to pass to initialization (the child positions).
799 * \param[in] scanner Scanner data structure.
800 * \returns Pointer to the created selection element.
802 * Implements _gmx_sel_init_keyword_evaluator() for \ref POS_VALUE input
805 static gmx::SelectionTreeElementPointer
806 init_evaluator_pos(gmx_ana_selmethod_t *method,
807 const gmx::SelectionParserParameterList ¶ms,
810 if ((method->flags & (SMETH_SINGLEVAL | SMETH_VARNUMVAL))
811 || method->outinit || method->pupdate == NULL)
814 = gmx::formatString("Keyword '%s' cannot be evaluated in this context",
816 GMX_THROW(gmx::InvalidInputError(message));
819 gmx::SelectionTreeElementPointer sel(
820 new gmx::SelectionTreeElement(
821 SEL_EXPRESSION, _gmx_sel_lexer_get_current_location(scanner)));
822 _gmx_selelem_set_method(sel, method, scanner);
824 t_methoddata_kweval *data
825 = new t_methoddata_kweval(sel->u.expr.method, sel->u.expr.mdata);
827 snew(sel->u.expr.method, 1);
828 sel->u.expr.method->name = data->kwmethod->name;
829 sel->u.expr.method->type = data->kwmethod->type;
830 sel->u.expr.method->flags = data->kwmethod->flags | SMETH_SINGLEVAL;
831 sel->u.expr.method->init_data = NULL;
832 sel->u.expr.method->set_poscoll = NULL;
833 sel->u.expr.method->init = method->init ? &init_kweval : NULL;
834 sel->u.expr.method->outinit = NULL;
835 sel->u.expr.method->free = &free_data_kweval;
836 sel->u.expr.method->init_frame = method->init_frame ? &init_frame_kweval : NULL;
837 sel->u.expr.method->update = &evaluate_kweval_pos;
838 sel->u.expr.method->pupdate = NULL;
839 sel->u.expr.method->nparams = asize(smparams_kweval_pos);
840 sel->u.expr.method->param = smparams_kweval_pos;
841 _gmx_selelem_init_method_params(sel, scanner);
842 sel->u.expr.mdata = data;
844 sel->u.expr.method->param[0].val.u.p = &data->p;
846 _gmx_sel_parse_params(params, sel->u.expr.method->nparams,
847 sel->u.expr.method->param, sel, scanner);
851 gmx::SelectionTreeElementPointer
852 _gmx_sel_init_keyword_evaluator(gmx_ana_selmethod_t *method,
853 const gmx::SelectionTreeElementPointer &child,
856 gmx::SelectionParserParameterList params;
858 gmx::SelectionParserParameter::createFromExpression(NULL, child));
859 if (child->v.type == GROUP_VALUE)
861 return init_evaluator_group(method, params, scanner);
863 else if (child->v.type == POS_VALUE)
865 return init_evaluator_pos(method, params, scanner);
869 std::string text(_gmx_sel_lexer_get_text(scanner, child->location()));
871 = gmx::formatString("Expression '%s' cannot be used to evaluate keywords",
873 GMX_THROW(gmx::InvalidInputError(message));