Additional index group checks for selections
[alexxy/gromacs.git] / src / gromacs / selection / selectioncollection.cpp
1 /*
2  * This file is part of the GROMACS molecular simulation package.
3  *
4  * Copyright (c) 2010,2011,2012,2013,2014, 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.
8  *
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.
13  *
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.
18  *
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.
23  *
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.
31  *
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.
34  */
35 /*! \internal \file
36  * \brief
37  * Implements gmx::SelectionCollection.
38  *
39  * \author Teemu Murtola <teemu.murtola@gmail.com>
40  * \ingroup module_selection
41  */
42 #include "selectioncollection.h"
43
44 #include <cctype>
45 #include <cstdio>
46
47 #include <string>
48 #include <vector>
49
50 #include <boost/shared_ptr.hpp>
51
52 #include "gromacs/legacyheaders/oenv.h"
53
54 #include "gromacs/onlinehelp/helpmanager.h"
55 #include "gromacs/onlinehelp/helpwritercontext.h"
56 #include "gromacs/options/basicoptions.h"
57 #include "gromacs/options/options.h"
58 #include "gromacs/selection/selection.h"
59 #include "gromacs/topology/topology.h"
60 #include "gromacs/utility/exceptions.h"
61 #include "gromacs/utility/file.h"
62 #include "gromacs/utility/gmxassert.h"
63 #include "gromacs/utility/messagestringcollector.h"
64 #include "gromacs/utility/smalloc.h"
65 #include "gromacs/utility/stringutil.h"
66
67 #include "compiler.h"
68 #include "mempool.h"
69 #include "parser.h"
70 #include "poscalc.h"
71 #include "scanner.h"
72 #include "selection.h"
73 #include "selectioncollection-impl.h"
74 #include "selelem.h"
75 #include "selhelp.h"
76 #include "selmethod.h"
77 #include "symrec.h"
78
79 namespace gmx
80 {
81
82 /********************************************************************
83  * SelectionCollection::Impl
84  */
85
86 SelectionCollection::Impl::Impl()
87     : debugLevel_(0), bExternalGroupsSet_(false), grps_(NULL)
88 {
89     sc_.nvars     = 0;
90     sc_.varstrs   = NULL;
91     sc_.top       = NULL;
92     gmx_ana_index_clear(&sc_.gall);
93     sc_.mempool   = NULL;
94     sc_.symtab.reset(new SelectionParserSymbolTable);
95     gmx_ana_selmethod_register_defaults(sc_.symtab.get());
96 }
97
98
99 SelectionCollection::Impl::~Impl()
100 {
101     clearSymbolTable();
102     // The tree must be freed before the SelectionData objects, since the
103     // tree may hold references to the position data in SelectionData.
104     sc_.root.reset();
105     sc_.sel.clear();
106     for (int i = 0; i < sc_.nvars; ++i)
107     {
108         sfree(sc_.varstrs[i]);
109     }
110     sfree(sc_.varstrs);
111     gmx_ana_index_deinit(&sc_.gall);
112     if (sc_.mempool)
113     {
114         _gmx_sel_mempool_destroy(sc_.mempool);
115     }
116 }
117
118
119 void
120 SelectionCollection::Impl::clearSymbolTable()
121 {
122     sc_.symtab.reset();
123 }
124
125
126 namespace
127 {
128
129 /*! \brief
130  * Reads a single selection line from stdin.
131  *
132  * \param[in]  infile        File to read from (typically File::standardInput()).
133  * \param[in]  bInteractive  Whether to print interactive prompts.
134  * \param[out] line          The read line in stored here.
135  * \returns true if something was read, false if at end of input.
136  *
137  * Handles line continuation, reading also the continuing line(s) in one call.
138  */
139 bool promptLine(File *infile, bool bInteractive, std::string *line)
140 {
141     if (bInteractive)
142     {
143         fprintf(stderr, "> ");
144     }
145     if (!infile->readLineWithTrailingSpace(line))
146     {
147         return false;
148     }
149     while (endsWith(*line, "\\\n"))
150     {
151         line->resize(line->length() - 2);
152         if (bInteractive)
153         {
154             fprintf(stderr, "... ");
155         }
156         std::string buffer;
157         // Return value ignored, buffer remains empty and works correctly
158         // if there is nothing to read.
159         infile->readLineWithTrailingSpace(&buffer);
160         line->append(buffer);
161     }
162     if (endsWith(*line, "\n"))
163     {
164         line->resize(line->length() - 1);
165     }
166     else if (bInteractive)
167     {
168         fprintf(stderr, "\n");
169     }
170     return true;
171 }
172
173 /*! \brief
174  * Helper function for tokenizing the input and pushing them to the parser.
175  *
176  * \param     scanner       Tokenizer data structure.
177  * \param     parserState   Parser data structure.
178  * \param[in] bInteractive  Whether to operate in interactive mode.
179  *
180  * Repeatedly reads tokens using \p scanner and pushes them to the parser with
181  * \p parserState until there is no more input, or until enough input is given
182  * (only in interactive mode).
183  */
184 int runParserLoop(yyscan_t scanner, _gmx_sel_yypstate *parserState,
185                   bool bInteractive)
186 {
187     int status    = YYPUSH_MORE;
188     int prevToken = 0;
189     do
190     {
191         YYSTYPE value;
192         int     token = _gmx_sel_yylex(&value, scanner);
193         if (bInteractive)
194         {
195             if (token == 0)
196             {
197                 break;
198             }
199             // Empty commands cause the interactive parser to print out
200             // status information. This avoids producing those unnecessarily,
201             // e.g., from "resname RA;;".
202             if (prevToken == CMD_SEP && token == CMD_SEP)
203             {
204                 continue;
205             }
206             prevToken = token;
207         }
208         status = _gmx_sel_yypush_parse(parserState, token, &value, scanner);
209     }
210     while (status == YYPUSH_MORE);
211     _gmx_sel_lexer_rethrow_exception_if_occurred(scanner);
212     return status;
213 }
214
215 /*! \brief
216  * Print current status in response to empty line in interactive input.
217  *
218  * \param[in] sc             Selection collection data structure.
219  * \param[in] grps           Available index groups.
220  * \param[in] firstSelection Index of first selection from this interactive
221  *     session.
222  * \param[in] maxCount       Maximum number of selections.
223  * \param[in] context        Context to print for what the selections are for.
224  * \param[in] bFirst         Whether this is the header that is printed before
225  *     any user input.
226  *
227  * Prints the available index groups and currently provided selections.
228  */
229 void printCurrentStatus(gmx_ana_selcollection_t *sc, gmx_ana_indexgrps_t *grps,
230                         size_t firstSelection, int maxCount,
231                         const std::string &context, bool bFirst)
232 {
233     if (grps != NULL)
234     {
235         std::fprintf(stderr, "Available static index groups:\n");
236         gmx_ana_indexgrps_print(stderr, grps, 0);
237     }
238     std::fprintf(stderr, "Specify ");
239     if (maxCount < 0)
240     {
241         std::fprintf(stderr, "any number of selections");
242     }
243     else if (maxCount == 1)
244     {
245         std::fprintf(stderr, "a selection");
246     }
247     else
248     {
249         std::fprintf(stderr, "%d selections", maxCount);
250     }
251     std::fprintf(stderr, "%s%s:\n",
252                  context.empty() ? "" : " ", context.c_str());
253     std::fprintf(stderr,
254                  "(one per line, <enter> for status/groups, 'help' for help%s)\n",
255                  maxCount < 0 ? ", Ctrl-D to end" : "");
256     if (!bFirst && (sc->nvars > 0 || sc->sel.size() > firstSelection))
257     {
258         std::fprintf(stderr, "Currently provided selections:\n");
259         for (int i = 0; i < sc->nvars; ++i)
260         {
261             std::fprintf(stderr, "     %s\n", sc->varstrs[i]);
262         }
263         for (size_t i = firstSelection; i < sc->sel.size(); ++i)
264         {
265             std::fprintf(stderr, " %2d. %s\n",
266                          static_cast<int>(i - firstSelection + 1),
267                          sc->sel[i]->selectionText());
268         }
269         if (maxCount > 0)
270         {
271             const int remaining
272                 = maxCount - static_cast<int>(sc->sel.size() - firstSelection);
273             std::fprintf(stderr, "(%d more selection%s required)\n",
274                          remaining, remaining > 1 ? "s" : "");
275         }
276     }
277 }
278
279 /*! \brief
280  * Prints selection help in interactive selection input.
281  *
282  * \param[in] sc    Selection collection data structure.
283  * \param[in] line  Line of user input requesting help (starting with `help`).
284  *
285  * Initializes the selection help if not yet initialized, and finds the help
286  * topic based on words on the input line.
287  */
288 void printHelp(gmx_ana_selcollection_t *sc, const std::string &line)
289 {
290     if (sc->rootHelp.get() == NULL)
291     {
292         sc->rootHelp = createSelectionHelpTopic();
293     }
294     HelpWriterContext context(&File::standardError(),
295                               eHelpOutputFormat_Console);
296     HelpManager       manager(*sc->rootHelp, context);
297     try
298     {
299         std::vector<std::string>                 topic = splitString(line);
300         std::vector<std::string>::const_iterator value;
301         // First item in the list is the 'help' token.
302         for (value = topic.begin() + 1; value != topic.end(); ++value)
303         {
304             manager.enterTopic(*value);
305         }
306     }
307     catch (const InvalidInputError &ex)
308     {
309         fprintf(stderr, "%s\n", ex.what());
310         return;
311     }
312     manager.writeCurrentTopic();
313 }
314
315 /*! \brief
316  * Helper function that runs the parser once the tokenizer has been
317  * initialized.
318  *
319  * \param[in,out] scanner Scanner data structure.
320  * \param[in]     bStdIn  Whether to use a line-based reading
321  *      algorithm designed for interactive input.
322  * \param[in]     maxnr   Maximum number of selections to parse
323  *      (if -1, parse as many as provided by the user).
324  * \param[in]     context Context to print for what the selections are for.
325  * \returns       Vector of parsed selections.
326  * \throws        std::bad_alloc if out of memory.
327  * \throws        InvalidInputError if there is a parsing error.
328  *
329  * Used internally to implement parseFromStdin(), parseFromFile() and
330  * parseFromString().
331  */
332 SelectionList runParser(yyscan_t scanner, bool bStdIn, int maxnr,
333                         const std::string &context)
334 {
335     boost::shared_ptr<void>  scannerGuard(scanner, &_gmx_sel_free_lexer);
336     gmx_ana_selcollection_t *sc   = _gmx_sel_lexer_selcollection(scanner);
337     gmx_ana_indexgrps_t     *grps = _gmx_sel_lexer_indexgrps(scanner);
338
339     MessageStringCollector   errors;
340     _gmx_sel_set_lexer_error_reporter(scanner, &errors);
341
342     size_t oldCount = sc->sel.size();
343     bool   bOk      = false;
344     {
345         boost::shared_ptr<_gmx_sel_yypstate> parserState(
346                 _gmx_sel_yypstate_new(), &_gmx_sel_yypstate_delete);
347         if (bStdIn)
348         {
349             File       &stdinFile(File::standardInput());
350             const bool  bInteractive = _gmx_sel_is_lexer_interactive(scanner);
351             if (bInteractive)
352             {
353                 printCurrentStatus(sc, grps, oldCount, maxnr, context, true);
354             }
355             std::string line;
356             int         status;
357             while (promptLine(&stdinFile, bInteractive, &line))
358             {
359                 if (bInteractive)
360                 {
361                     line = stripString(line);
362                     if (line.empty())
363                     {
364                         printCurrentStatus(sc, grps, oldCount, maxnr, context, false);
365                         continue;
366                     }
367                     if (startsWith(line, "help")
368                         && (line[4] == 0 || std::isspace(line[4])))
369                     {
370                         printHelp(sc, line);
371                         continue;
372                     }
373                 }
374                 line.append("\n");
375                 _gmx_sel_set_lex_input_str(scanner, line.c_str());
376                 status = runParserLoop(scanner, parserState.get(), true);
377                 if (status != YYPUSH_MORE)
378                 {
379                     // TODO: Check if there is more input, and issue an
380                     // error/warning if some input was ignored.
381                     goto early_termination;
382                 }
383                 if (!errors.isEmpty() && bInteractive)
384                 {
385                     fprintf(stderr, "%s", errors.toString().c_str());
386                     errors.clear();
387                 }
388             }
389             status = _gmx_sel_yypush_parse(parserState.get(), 0, NULL,
390                                            scanner);
391             _gmx_sel_lexer_rethrow_exception_if_occurred(scanner);
392 early_termination:
393             bOk = (status == 0);
394         }
395         else
396         {
397             int status = runParserLoop(scanner, parserState.get(), false);
398             bOk = (status == 0);
399         }
400     }
401     scannerGuard.reset();
402     int nr = sc->sel.size() - oldCount;
403     if (maxnr > 0 && nr != maxnr)
404     {
405         bOk = false;
406         errors.append("Too few selections provided");
407     }
408
409     // TODO: Remove added selections from the collection if parsing failed?
410     if (!bOk || !errors.isEmpty())
411     {
412         GMX_ASSERT(!bOk && !errors.isEmpty(), "Inconsistent error reporting");
413         GMX_THROW(InvalidInputError(errors.toString()));
414     }
415
416     SelectionList                     result;
417     SelectionDataList::const_iterator i;
418     result.reserve(nr);
419     for (i = sc->sel.begin() + oldCount; i != sc->sel.end(); ++i)
420     {
421         result.push_back(Selection(i->get()));
422     }
423     return result;
424 }
425
426 /*! \brief
427  * Checks that index groups have valid atom indices.
428  *
429  * \param[in]    root    Root of selection tree to process.
430  * \param[in]    natoms  Maximum number of atoms that the selections are set
431  *     to evaluate.
432  * \param        errors  Object for reporting any error messages.
433  * \throws std::bad_alloc if out of memory.
434  *
435  * Recursively checks the selection tree for index groups.
436  * Each found group is checked that it only contains atom indices that match
437  * the topology/maximum number of atoms set for the selection collection.
438  * Any issues are reported to \p errors.
439  */
440 void checkExternalGroups(const SelectionTreeElementPointer &root,
441                          int                                natoms,
442                          ExceptionInitializer              *errors)
443 {
444     if (root->type == SEL_CONST && root->v.type == GROUP_VALUE)
445     {
446         try
447         {
448             root->checkIndexGroup(natoms);
449         }
450         catch (const UserInputError &)
451         {
452             errors->addCurrentExceptionAsNested();
453         }
454     }
455
456     SelectionTreeElementPointer child = root->child;
457     while (child)
458     {
459         checkExternalGroups(child, natoms, errors);
460         child = child->next;
461     }
462 }
463
464 }   // namespace
465
466
467 void SelectionCollection::Impl::resolveExternalGroups(
468         const SelectionTreeElementPointer &root,
469         ExceptionInitializer              *errors)
470 {
471
472     if (root->type == SEL_GROUPREF)
473     {
474         try
475         {
476             root->resolveIndexGroupReference(grps_, sc_.gall.isize);
477         }
478         catch (const UserInputError &)
479         {
480             errors->addCurrentExceptionAsNested();
481         }
482     }
483
484     SelectionTreeElementPointer child = root->child;
485     while (child)
486     {
487         resolveExternalGroups(child, errors);
488         root->flags |= (child->flags & SEL_UNSORTED);
489         child        = child->next;
490     }
491 }
492
493
494 /********************************************************************
495  * SelectionCollection
496  */
497
498 SelectionCollection::SelectionCollection()
499     : impl_(new Impl)
500 {
501 }
502
503
504 SelectionCollection::~SelectionCollection()
505 {
506 }
507
508
509 void
510 SelectionCollection::initOptions(Options *options)
511 {
512     const char * const debug_levels[]
513         = { "no", "basic", "compile", "eval", "full" };
514
515     bool bAllowNonAtomOutput = false;
516     SelectionDataList::const_iterator iter;
517     for (iter = impl_->sc_.sel.begin(); iter != impl_->sc_.sel.end(); ++iter)
518     {
519         const internal::SelectionData &sel = **iter;
520         if (!sel.hasFlag(efSelection_OnlyAtoms))
521         {
522             bAllowNonAtomOutput = true;
523         }
524     }
525
526     const char *const *postypes = PositionCalculationCollection::typeEnumValues;
527     options->addOption(StringOption("selrpos")
528                            .enumValueFromNullTerminatedArray(postypes)
529                            .store(&impl_->rpost_).defaultValue(postypes[0])
530                            .description("Selection reference positions"));
531     if (bAllowNonAtomOutput)
532     {
533         options->addOption(StringOption("seltype")
534                                .enumValueFromNullTerminatedArray(postypes)
535                                .store(&impl_->spost_).defaultValue(postypes[0])
536                                .description("Default selection output positions"));
537     }
538     else
539     {
540         impl_->spost_ = postypes[0];
541     }
542     GMX_RELEASE_ASSERT(impl_->debugLevel_ >= 0 && impl_->debugLevel_ <= 4,
543                        "Debug level out of range");
544     options->addOption(StringOption("seldebug").hidden(impl_->debugLevel_ == 0)
545                            .enumValue(debug_levels)
546                            .defaultValue(debug_levels[impl_->debugLevel_])
547                            .storeEnumIndex(&impl_->debugLevel_)
548                            .description("Print out selection trees for debugging"));
549 }
550
551
552 void
553 SelectionCollection::setReferencePosType(const char *type)
554 {
555     GMX_RELEASE_ASSERT(type != NULL, "Cannot assign NULL position type");
556     // Check that the type is valid, throw if it is not.
557     e_poscalc_t  dummytype;
558     int          dummyflags;
559     PositionCalculationCollection::typeFromEnum(type, &dummytype, &dummyflags);
560     impl_->rpost_ = type;
561 }
562
563
564 void
565 SelectionCollection::setOutputPosType(const char *type)
566 {
567     GMX_RELEASE_ASSERT(type != NULL, "Cannot assign NULL position type");
568     // Check that the type is valid, throw if it is not.
569     e_poscalc_t  dummytype;
570     int          dummyflags;
571     PositionCalculationCollection::typeFromEnum(type, &dummytype, &dummyflags);
572     impl_->spost_ = type;
573 }
574
575
576 void
577 SelectionCollection::setDebugLevel(int debugLevel)
578 {
579     impl_->debugLevel_ = debugLevel;
580 }
581
582
583 void
584 SelectionCollection::setTopology(t_topology *top, int natoms)
585 {
586     GMX_RELEASE_ASSERT(natoms > 0 || top != NULL,
587                        "The number of atoms must be given if there is no topology");
588     // Get the number of atoms from the topology if it is not given.
589     if (natoms <= 0)
590     {
591         natoms = top->atoms.nr;
592     }
593     if (impl_->bExternalGroupsSet_)
594     {
595         ExceptionInitializer        errors("Invalid index group references encountered");
596         SelectionTreeElementPointer root = impl_->sc_.root;
597         while (root)
598         {
599             checkExternalGroups(root, natoms, &errors);
600             root = root->next;
601         }
602         if (errors.hasNestedExceptions())
603         {
604             GMX_THROW(InconsistentInputError(errors));
605         }
606     }
607     gmx_ana_selcollection_t *sc = &impl_->sc_;
608     // Do this first, as it allocates memory, while the others don't throw.
609     gmx_ana_index_init_simple(&sc->gall, natoms);
610     sc->pcc.setTopology(top);
611     sc->top = top;
612 }
613
614
615 void
616 SelectionCollection::setIndexGroups(gmx_ana_indexgrps_t *grps)
617 {
618     GMX_RELEASE_ASSERT(grps == NULL || !impl_->bExternalGroupsSet_,
619                        "Can only set external groups once or clear them afterwards");
620     impl_->grps_               = grps;
621     impl_->bExternalGroupsSet_ = true;
622
623     ExceptionInitializer        errors("Invalid index group reference(s)");
624     SelectionTreeElementPointer root = impl_->sc_.root;
625     while (root)
626     {
627         impl_->resolveExternalGroups(root, &errors);
628         root->checkUnsortedAtoms(true, &errors);
629         root = root->next;
630     }
631     if (errors.hasNestedExceptions())
632     {
633         GMX_THROW(InconsistentInputError(errors));
634     }
635     for (size_t i = 0; i < impl_->sc_.sel.size(); ++i)
636     {
637         impl_->sc_.sel[i]->refreshName();
638     }
639 }
640
641
642 bool
643 SelectionCollection::requiresTopology() const
644 {
645     e_poscalc_t  type;
646     int          flags;
647
648     if (!impl_->rpost_.empty())
649     {
650         flags = 0;
651         // Should not throw, because has been checked earlier.
652         PositionCalculationCollection::typeFromEnum(impl_->rpost_.c_str(),
653                                                     &type, &flags);
654         if (type != POS_ATOM)
655         {
656             return true;
657         }
658     }
659     if (!impl_->spost_.empty())
660     {
661         flags = 0;
662         // Should not throw, because has been checked earlier.
663         PositionCalculationCollection::typeFromEnum(impl_->spost_.c_str(),
664                                                     &type, &flags);
665         if (type != POS_ATOM)
666         {
667             return true;
668         }
669     }
670
671     SelectionTreeElementPointer sel = impl_->sc_.root;
672     while (sel)
673     {
674         if (_gmx_selelem_requires_top(*sel))
675         {
676             return true;
677         }
678         sel = sel->next;
679     }
680     return false;
681 }
682
683
684 SelectionList
685 SelectionCollection::parseFromStdin(int nr, bool bInteractive,
686                                     const std::string &context)
687 {
688     yyscan_t scanner;
689
690     _gmx_sel_init_lexer(&scanner, &impl_->sc_, bInteractive, nr,
691                         impl_->bExternalGroupsSet_,
692                         impl_->grps_);
693     return runParser(scanner, true, nr, context);
694 }
695
696
697 SelectionList
698 SelectionCollection::parseFromFile(const std::string &filename)
699 {
700
701     try
702     {
703         yyscan_t scanner;
704         File     file(filename, "r");
705         // TODO: Exception-safe way of using the lexer.
706         _gmx_sel_init_lexer(&scanner, &impl_->sc_, false, -1,
707                             impl_->bExternalGroupsSet_,
708                             impl_->grps_);
709         _gmx_sel_set_lex_input_file(scanner, file.handle());
710         return runParser(scanner, false, -1, std::string());
711     }
712     catch (GromacsException &ex)
713     {
714         ex.prependContext(formatString(
715                                   "Error in parsing selections from file '%s'",
716                                   filename.c_str()));
717         throw;
718     }
719 }
720
721
722 SelectionList
723 SelectionCollection::parseFromString(const std::string &str)
724 {
725     yyscan_t scanner;
726
727     _gmx_sel_init_lexer(&scanner, &impl_->sc_, false, -1,
728                         impl_->bExternalGroupsSet_,
729                         impl_->grps_);
730     _gmx_sel_set_lex_input_str(scanner, str.c_str());
731     return runParser(scanner, false, -1, std::string());
732 }
733
734
735 void
736 SelectionCollection::compile()
737 {
738     if (impl_->sc_.top == NULL && requiresTopology())
739     {
740         GMX_THROW(InconsistentInputError("Selection requires topology information, but none provided"));
741     }
742     if (!impl_->bExternalGroupsSet_)
743     {
744         setIndexGroups(NULL);
745     }
746     if (impl_->debugLevel_ >= 1)
747     {
748         printTree(stderr, false);
749     }
750
751     SelectionCompiler compiler;
752     compiler.compile(this);
753
754     if (impl_->debugLevel_ >= 1)
755     {
756         std::fprintf(stderr, "\n");
757         printTree(stderr, false);
758         std::fprintf(stderr, "\n");
759         impl_->sc_.pcc.printTree(stderr);
760         std::fprintf(stderr, "\n");
761     }
762     impl_->sc_.pcc.initEvaluation();
763     if (impl_->debugLevel_ >= 1)
764     {
765         impl_->sc_.pcc.printTree(stderr);
766         std::fprintf(stderr, "\n");
767     }
768
769     // TODO: It would be nicer to associate the name of the selection option
770     // (if available) to the error message.
771     SelectionDataList::const_iterator iter;
772     for (iter = impl_->sc_.sel.begin(); iter != impl_->sc_.sel.end(); ++iter)
773     {
774         const internal::SelectionData &sel = **iter;
775         if (sel.hasFlag(efSelection_OnlyAtoms))
776         {
777             if (sel.type() != INDEX_ATOM)
778             {
779                 std::string message = formatString(
780                             "Selection '%s' does not evaluate to individual atoms. "
781                             "This is not allowed in this context.",
782                             sel.selectionText());
783                 GMX_THROW(InvalidInputError(message));
784             }
785         }
786         if (sel.hasFlag(efSelection_DisallowEmpty))
787         {
788             if (sel.posCount() == 0)
789             {
790                 std::string message = formatString(
791                             "Selection '%s' never matches any atoms.",
792                             sel.selectionText());
793                 GMX_THROW(InvalidInputError(message));
794             }
795         }
796     }
797 }
798
799
800 void
801 SelectionCollection::evaluate(t_trxframe *fr, t_pbc *pbc)
802 {
803     impl_->sc_.pcc.initFrame();
804
805     SelectionEvaluator evaluator;
806     evaluator.evaluate(this, fr, pbc);
807
808     if (impl_->debugLevel_ >= 3)
809     {
810         std::fprintf(stderr, "\n");
811         printTree(stderr, true);
812     }
813 }
814
815
816 void
817 SelectionCollection::evaluateFinal(int nframes)
818 {
819     SelectionEvaluator evaluator;
820     evaluator.evaluateFinal(this, nframes);
821 }
822
823
824 void
825 SelectionCollection::printTree(FILE *fp, bool bValues) const
826 {
827     SelectionTreeElementPointer sel = impl_->sc_.root;
828     while (sel)
829     {
830         _gmx_selelem_print_tree(fp, *sel, bValues, 0);
831         sel = sel->next;
832     }
833 }
834
835
836 void
837 SelectionCollection::printXvgrInfo(FILE *out, output_env_t oenv) const
838 {
839     if (output_env_get_xvg_format(oenv) != exvgNONE)
840     {
841         const gmx_ana_selcollection_t &sc = impl_->sc_;
842         std::fprintf(out, "# Selections:\n");
843         for (int i = 0; i < sc.nvars; ++i)
844         {
845             std::fprintf(out, "#   %s\n", sc.varstrs[i]);
846         }
847         for (size_t i = 0; i < sc.sel.size(); ++i)
848         {
849             std::fprintf(out, "#   %s\n", sc.sel[i]->selectionText());
850         }
851         std::fprintf(out, "#\n");
852     }
853 }
854
855 } // namespace gmx