src/gromacs/selection/selectioncollection.cpp

   1 /*
   2  *
   3  *                This source code is part of
   4  *
   5  *                 G   R   O   M   A   C   S
   6  *
   7  *          GROningen MAchine for Chemical Simulations
   8  *
   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.
  13
  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.
  18  *
  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.
  25  *
  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.
  28  *
  29  * For more info, check our website at http://www.gromacs.org
  30  */
  31 /*! \internal \file
  32  * \brief
  33  * Implements gmx::SelectionCollection.
  34  *
  35  * \author Teemu Murtola <teemu.murtola@cbr.su.se>
  36  * \ingroup module_selection
  37  */
  38 #ifdef HAVE_CONFIG_H
  39 #include <config.h>
  40 #endif
  41
  42 #include <cstdio>
  43
  44 #include <boost/shared_ptr.hpp>
  45
  46 #include "gromacs/legacyheaders/oenv.h"
  47 #include "gromacs/legacyheaders/smalloc.h"
  48 #include "gromacs/legacyheaders/xvgr.h"
  49
  50 #include "gromacs/options/basicoptions.h"
  51 #include "gromacs/options/options.h"
  52 #include "gromacs/selection/selection.h"
  53 #include "gromacs/selection/selectioncollection.h"
  54 #include "gromacs/utility/exceptions.h"
  55 #include "gromacs/utility/file.h"
  56 #include "gromacs/utility/gmxassert.h"
  57 #include "gromacs/utility/messagestringcollector.h"
  58 #include "gromacs/utility/stringutil.h"
  59
  60 #include "compiler.h"
  61 #include "mempool.h"
  62 #include "parser.h"
  63 #include "poscalc.h"
  64 #include "scanner.h"
  65 #include "selectioncollection-impl.h"
  66 #include "selelem.h"
  67 #include "selhelp.h"
  68 #include "selmethod.h"
  69 #include "symrec.h"
  70
  71 namespace gmx
  72 {
  73
  74 /********************************************************************
  75  * SelectionCollection::Impl
  76  */
  77
  78 SelectionCollection::Impl::Impl()
  79     : debugLevel_(0), bExternalGroupsSet_(false), grps_(NULL)
  80 {
  81     sc_.nvars     = 0;
  82     sc_.varstrs   = NULL;
  83     sc_.top       = NULL;
  84     gmx_ana_index_clear(&sc_.gall);
  85     sc_.mempool   = NULL;
  86     sc_.symtab.reset(new SelectionParserSymbolTable);
  87     gmx_ana_selmethod_register_defaults(sc_.symtab.get());
  88 }
  89
  90
  91 SelectionCollection::Impl::~Impl()
  92 {
  93     clearSymbolTable();
  94     // The tree must be freed before the SelectionData objects, since the
  95     // tree may hold references to the position data in SelectionData.
  96     sc_.root.reset();
  97     sc_.sel.clear();
  98     for (int i = 0; i < sc_.nvars; ++i)
  99     {
 100         sfree(sc_.varstrs[i]);
 101     }
 102     sfree(sc_.varstrs);
 103     gmx_ana_index_deinit(&sc_.gall);
 104     if (sc_.mempool)
 105     {
 106         _gmx_sel_mempool_destroy(sc_.mempool);
 107     }
 108 }
 109
 110
 111 void
 112 SelectionCollection::Impl::clearSymbolTable()
 113 {
 114     sc_.symtab.reset();
 115 }
 116
 117
 118 namespace
 119 {
 120
 121 /*! \brief
 122  * Reads a single selection line from stdin.
 123  *
 124  * \param[in]  infile        File to read from (typically File::standardInput()).
 125  * \param[in]  bInteractive  Whether to print interactive prompts.
 126  * \param[out] line          The read line in stored here.
 127  * \returns true if something was read, false if at end of input.
 128  *
 129  * Handles line continuation, reading also the continuing line(s) in one call.
 130  */
 131 bool promptLine(File *infile, bool bInteractive, std::string *line)
 132 {
 133     if (bInteractive)
 134     {
 135         fprintf(stderr, "> ");
 136     }
 137     if (!infile->readLine(line))
 138     {
 139         return false;
 140     }
 141     while(endsWith(*line, "\\\n"))
 142     {
 143         line->resize(line->length() - 2);
 144         if (bInteractive)
 145         {
 146             fprintf(stderr, "... ");
 147         }
 148         std::string buffer;
 149         // Return value ignored, buffer remains empty and works correctly
 150         // if there is nothing to read.
 151         infile->readLine(&buffer);
 152         line->append(buffer);
 153     }
 154     if (endsWith(*line, "\n"))
 155     {
 156         line->resize(line->length() - 1);
 157     }
 158     else if (bInteractive)
 159     {
 160         fprintf(stderr, "\n");
 161     }
 162     return true;
 163 }
 164
 165 /*! \brief
 166  * Helper function for tokenizing the input and pushing them to the parser.
 167  *
 168  * \param     scanner       Tokenizer data structure.
 169  * \param     parserState   Parser data structure.
 170  * \param[in] bInteractive  Whether to operate in interactive mode.
 171  *
 172  * Repeatedly reads tokens using \p scanner and pushes them to the parser with
 173  * \p parserState until there is no more input, or until enough input is given
 174  * (only in interactive mode).
 175  */
 176 int runParserLoop(yyscan_t scanner, _gmx_sel_yypstate *parserState,
 177                   bool bInteractive)
 178 {
 179     int status = YYPUSH_MORE;
 180     int prevToken = 0;
 181     do
 182     {
 183         YYSTYPE value;
 184         int token = _gmx_sel_yylex(&value, scanner);
 185         if (bInteractive)
 186         {
 187             if (token == 0)
 188             {
 189                 break;
 190             }
 191             // Empty commands cause the interactive parser to print out
 192             // status information. This avoids producing those unnecessarily,
 193             // e.g., from "resname RA;;".
 194             if (prevToken == CMD_SEP && token == CMD_SEP)
 195             {
 196                 continue;
 197             }
 198             prevToken = token;
 199         }
 200         status = _gmx_sel_yypush_parse(parserState, token, &value, scanner);
 201     }
 202     while (status == YYPUSH_MORE);
 203     _gmx_sel_lexer_rethrow_exception_if_occurred(scanner);
 204     return status;
 205 }
 206
 207 /*! \brief
 208  * Helper function that runs the parser once the tokenizer has been
 209  * initialized.
 210  *
 211  * \param[in,out] scanner Scanner data structure.
 212  * \param[in]     bStdIn  Whether to use a line-based reading
 213  *      algorithm designed for interactive input.
 214  * \param[in]     maxnr   Maximum number of selections to parse
 215  *      (if -1, parse as many as provided by the user).
 216  * \returns       Vector of parsed selections.
 217  * \throws        std::bad_alloc if out of memory.
 218  * \throws        InvalidInputError if there is a parsing error.
 219  *
 220  * Used internally to implement parseFromStdin(), parseFromFile() and
 221  * parseFromString().
 222  */
 223 SelectionList runParser(yyscan_t scanner, bool bStdIn, int maxnr)
 224 {
 225     boost::shared_ptr<void> scannerGuard(scanner, &_gmx_sel_free_lexer);
 226     gmx_ana_selcollection_t *sc = _gmx_sel_lexer_selcollection(scanner);
 227
 228     MessageStringCollector errors;
 229     _gmx_sel_set_lexer_error_reporter(scanner, &errors);
 230
 231     int oldCount = sc->sel.size();
 232     bool bOk = false;
 233     {
 234         boost::shared_ptr<_gmx_sel_yypstate> parserState(
 235                 _gmx_sel_yypstate_new(), &_gmx_sel_yypstate_delete);
 236         if (bStdIn)
 237         {
 238             File &stdinFile(File::standardInput());
 239             bool bInteractive = _gmx_sel_is_lexer_interactive(scanner);
 240             std::string line;
 241             int status;
 242             while (promptLine(&stdinFile, bInteractive, &line))
 243             {
 244                 line.append("\n");
 245                 _gmx_sel_set_lex_input_str(scanner, line.c_str());
 246                 status = runParserLoop(scanner, parserState.get(), true);
 247                 if (status != YYPUSH_MORE)
 248                 {
 249                     // TODO: Check if there is more input, and issue an
 250                     // error/warning if some input was ignored.
 251                     goto early_termination;
 252                 }
 253                 if (!errors.isEmpty() && bInteractive)
 254                 {
 255                     fprintf(stderr, "%s", errors.toString().c_str());
 256                     errors.clear();
 257                 }
 258             }
 259             status = _gmx_sel_yypush_parse(parserState.get(), 0, NULL,
 260                                            scanner);
 261             _gmx_sel_lexer_rethrow_exception_if_occurred(scanner);
 262 early_termination:
 263             bOk = (status == 0);
 264         }
 265         else
 266         {
 267             int status = runParserLoop(scanner, parserState.get(), false);
 268             bOk = (status == 0);
 269         }
 270     }
 271     scannerGuard.reset();
 272     int nr = sc->sel.size() - oldCount;
 273     if (maxnr > 0 && nr != maxnr)
 274     {
 275         bOk = false;
 276         errors.append("Too few selections provided");
 277     }
 278
 279     // TODO: Remove added selections from the collection if parsing failed?
 280     if (!bOk || !errors.isEmpty())
 281     {
 282         GMX_ASSERT(!bOk && !errors.isEmpty(), "Inconsistent error reporting");
 283         GMX_THROW(InvalidInputError(errors.toString()));
 284     }
 285
 286     SelectionList result;
 287     SelectionDataList::const_iterator i;
 288     result.reserve(nr);
 289     for (i = sc->sel.begin() + oldCount; i != sc->sel.end(); ++i)
 290     {
 291         result.push_back(Selection(i->get()));
 292     }
 293     return result;
 294 }
 295
 296 } // namespace
 297
 298
 299 void SelectionCollection::Impl::resolveExternalGroups(
 300         const SelectionTreeElementPointer &root,
 301         MessageStringCollector *errors)
 302 {
 303
 304     if (root->type == SEL_GROUPREF)
 305     {
 306         bool bOk = true;
 307         if (grps_ == NULL)
 308         {
 309             // TODO: Improve error messages
 310             errors->append("Unknown group referenced in a selection");
 311             bOk = false;
 312         }
 313         else if (root->u.gref.name != NULL)
 314         {
 315             char *name = root->u.gref.name;
 316             if (!gmx_ana_indexgrps_find(&root->u.cgrp, grps_, name))
 317             {
 318                 // TODO: Improve error messages
 319                 errors->append("Unknown group referenced in a selection");
 320                 bOk = false;
 321             }
 322             else
 323             {
 324                 sfree(name);
 325             }
 326         }
 327         else
 328         {
 329             if (!gmx_ana_indexgrps_extract(&root->u.cgrp, grps_,
 330                                            root->u.gref.id))
 331             {
 332                 // TODO: Improve error messages
 333                 errors->append("Unknown group referenced in a selection");
 334                 bOk = false;
 335             }
 336         }
 337         if (bOk)
 338         {
 339             root->type = SEL_CONST;
 340             root->setName(root->u.cgrp.name);
 341         }
 342     }
 343
 344     SelectionTreeElementPointer child = root->child;
 345     while (child)
 346     {
 347         resolveExternalGroups(child, errors);
 348         child = child->next;
 349     }
 350 }
 351
 352
 353 /********************************************************************
 354  * SelectionCollection
 355  */
 356
 357 SelectionCollection::SelectionCollection()
 358     : impl_(new Impl)
 359 {
 360 }
 361
 362
 363 SelectionCollection::~SelectionCollection()
 364 {
 365 }
 366
 367
 368 void
 369 SelectionCollection::initOptions(Options *options)
 370 {
 371     static const char * const debug_levels[]
 372         = {"no", "basic", "compile", "eval", "full", NULL};
 373     /*
 374     static const char * const desc[] = {
 375         "This program supports selections in addition to traditional",
 376         "index files. Use [TT]-select help[tt] for additional information,",
 377         "or type 'help' in the selection prompt.",
 378         NULL,
 379     };
 380     options.setDescription(desc);
 381     */
 382
 383     const char *const *postypes = PositionCalculationCollection::typeEnumValues;
 384     options->addOption(StringOption("selrpos").enumValue(postypes)
 385                            .store(&impl_->rpost_).defaultValue(postypes[0])
 386                            .description("Selection reference positions"));
 387     options->addOption(StringOption("seltype").enumValue(postypes)
 388                            .store(&impl_->spost_).defaultValue(postypes[0])
 389                            .description("Default selection output positions"));
 390     GMX_RELEASE_ASSERT(impl_->debugLevel_ >= 0 && impl_->debugLevel_ <= 4,
 391                        "Debug level out of range");
 392     options->addOption(StringOption("seldebug").hidden(impl_->debugLevel_ == 0)
 393                            .enumValue(debug_levels)
 394                            .defaultValue(debug_levels[impl_->debugLevel_])
 395                            .storeEnumIndex(&impl_->debugLevel_)
 396                            .description("Print out selection trees for debugging"));
 397 }
 398
 399
 400 void
 401 SelectionCollection::setReferencePosType(const char *type)
 402 {
 403     GMX_RELEASE_ASSERT(type != NULL, "Cannot assign NULL position type");
 404     // Check that the type is valid, throw if it is not.
 405     e_poscalc_t  dummytype;
 406     int          dummyflags;
 407     PositionCalculationCollection::typeFromEnum(type, &dummytype, &dummyflags);
 408     impl_->rpost_ = type;
 409 }
 410
 411
 412 void
 413 SelectionCollection::setOutputPosType(const char *type)
 414 {
 415     GMX_RELEASE_ASSERT(type != NULL, "Cannot assign NULL position type");
 416     // Check that the type is valid, throw if it is not.
 417     e_poscalc_t  dummytype;
 418     int          dummyflags;
 419     PositionCalculationCollection::typeFromEnum(type, &dummytype, &dummyflags);
 420     impl_->spost_ = type;
 421 }
 422
 423
 424 void
 425 SelectionCollection::setDebugLevel(int debugLevel)
 426 {
 427     impl_->debugLevel_ = debugLevel;
 428 }
 429
 430
 431 void
 432 SelectionCollection::setTopology(t_topology *top, int natoms)
 433 {
 434     GMX_RELEASE_ASSERT(natoms > 0 || top != NULL,
 435         "The number of atoms must be given if there is no topology");
 436     // Get the number of atoms from the topology if it is not given.
 437     if (natoms <= 0)
 438     {
 439         natoms = top->atoms.nr;
 440     }
 441     gmx_ana_selcollection_t *sc = &impl_->sc_;
 442     // Do this first, as it allocates memory, while the others don't throw.
 443     gmx_ana_index_init_simple(&sc->gall, natoms, NULL);
 444     sc->pcc.setTopology(top);
 445     sc->top = top;
 446 }
 447
 448
 449 void
 450 SelectionCollection::setIndexGroups(gmx_ana_indexgrps_t *grps)
 451 {
 452     GMX_RELEASE_ASSERT(grps == NULL || !impl_->bExternalGroupsSet_,
 453                        "Can only set external groups once or clear them afterwards");
 454     impl_->grps_ = grps;
 455     impl_->bExternalGroupsSet_ = true;
 456
 457     MessageStringCollector errors;
 458     SelectionTreeElementPointer root = impl_->sc_.root;
 459     while (root)
 460     {
 461         impl_->resolveExternalGroups(root, &errors);
 462         root = root->next;
 463     }
 464     if (!errors.isEmpty())
 465     {
 466         GMX_THROW(InvalidInputError(errors.toString()));
 467     }
 468 }
 469
 470
 471 bool
 472 SelectionCollection::requiresTopology() const
 473 {
 474     e_poscalc_t  type;
 475     int          flags;
 476
 477     if (!impl_->rpost_.empty())
 478     {
 479         flags = 0;
 480         // Should not throw, because has been checked earlier.
 481         PositionCalculationCollection::typeFromEnum(impl_->rpost_.c_str(),
 482                                                     &type, &flags);
 483         if (type != POS_ATOM)
 484         {
 485             return true;
 486         }
 487     }
 488     if (!impl_->spost_.empty())
 489     {
 490         flags = 0;
 491         // Should not throw, because has been checked earlier.
 492         PositionCalculationCollection::typeFromEnum(impl_->spost_.c_str(),
 493                                                     &type, &flags);
 494         if (type != POS_ATOM)
 495         {
 496             return true;
 497         }
 498     }
 499
 500     SelectionTreeElementPointer sel = impl_->sc_.root;
 501     while (sel)
 502     {
 503         if (_gmx_selelem_requires_top(*sel))
 504         {
 505             return true;
 506         }
 507         sel = sel->next;
 508     }
 509     return false;
 510 }
 511
 512
 513 SelectionList
 514 SelectionCollection::parseFromStdin(int nr, bool bInteractive)
 515 {
 516     yyscan_t scanner;
 517
 518     _gmx_sel_init_lexer(&scanner, &impl_->sc_, bInteractive, nr,
 519                         impl_->bExternalGroupsSet_,
 520                         impl_->grps_);
 521     return runParser(scanner, true, nr);
 522 }
 523
 524
 525 SelectionList
 526 SelectionCollection::parseFromFile(const std::string &filename)
 527 {
 528
 529     try
 530     {
 531         yyscan_t scanner;
 532         File file(filename, "r");
 533         // TODO: Exception-safe way of using the lexer.
 534         _gmx_sel_init_lexer(&scanner, &impl_->sc_, false, -1,
 535                             impl_->bExternalGroupsSet_,
 536                             impl_->grps_);
 537         _gmx_sel_set_lex_input_file(scanner, file.handle());
 538         return runParser(scanner, false, -1);
 539     }
 540     catch (GromacsException &ex)
 541     {
 542         ex.prependContext(formatString(
 543                     "Error in parsing selections from file '%s'",
 544                     filename.c_str()));
 545         throw;
 546     }
 547 }
 548
 549
 550 SelectionList
 551 SelectionCollection::parseFromString(const std::string &str)
 552 {
 553     yyscan_t scanner;
 554
 555     _gmx_sel_init_lexer(&scanner, &impl_->sc_, false, -1,
 556                         impl_->bExternalGroupsSet_,
 557                         impl_->grps_);
 558     _gmx_sel_set_lex_input_str(scanner, str.c_str());
 559     return runParser(scanner, false, -1);
 560 }
 561
 562
 563 void
 564 SelectionCollection::compile()
 565 {
 566     if (impl_->sc_.top == NULL && requiresTopology())
 567     {
 568         GMX_THROW(InconsistentInputError("Selection requires topology information, but none provided"));
 569     }
 570     if (!impl_->bExternalGroupsSet_)
 571     {
 572         setIndexGroups(NULL);
 573     }
 574     if (impl_->debugLevel_ >= 1)
 575     {
 576         printTree(stderr, false);
 577     }
 578
 579     SelectionCompiler compiler;
 580     compiler.compile(this);
 581
 582     if (impl_->debugLevel_ >= 1)
 583     {
 584         std::fprintf(stderr, "\n");
 585         printTree(stderr, false);
 586         std::fprintf(stderr, "\n");
 587         impl_->sc_.pcc.printTree(stderr);
 588         std::fprintf(stderr, "\n");
 589     }
 590     impl_->sc_.pcc.initEvaluation();
 591     if (impl_->debugLevel_ >= 1)
 592     {
 593         impl_->sc_.pcc.printTree(stderr);
 594         std::fprintf(stderr, "\n");
 595     }
 596 }
 597
 598
 599 void
 600 SelectionCollection::evaluate(t_trxframe *fr, t_pbc *pbc)
 601 {
 602     impl_->sc_.pcc.initFrame();
 603
 604     SelectionEvaluator evaluator;
 605     evaluator.evaluate(this, fr, pbc);
 606
 607     if (impl_->debugLevel_ >= 3)
 608     {
 609         std::fprintf(stderr, "\n");
 610         printTree(stderr, true);
 611     }
 612 }
 613
 614
 615 void
 616 SelectionCollection::evaluateFinal(int nframes)
 617 {
 618     SelectionEvaluator evaluator;
 619     evaluator.evaluateFinal(this, nframes);
 620 }
 621
 622
 623 void
 624 SelectionCollection::printTree(FILE *fp, bool bValues) const
 625 {
 626     SelectionTreeElementPointer sel = impl_->sc_.root;
 627     while (sel)
 628     {
 629         _gmx_selelem_print_tree(fp, *sel, bValues, 0);
 630         sel = sel->next;
 631     }
 632 }
 633
 634
 635 void
 636 SelectionCollection::printXvgrInfo(FILE *out, output_env_t oenv) const
 637 {
 638     if (output_env_get_xvg_format(oenv) != exvgNONE)
 639     {
 640         const gmx_ana_selcollection_t &sc = impl_->sc_;
 641         std::fprintf(out, "# Selections:\n");
 642         for (int i = 0; i < sc.nvars; ++i)
 643         {
 644             std::fprintf(out, "#   %s\n", sc.varstrs[i]);
 645         }
 646         for (size_t i = 0; i < sc.sel.size(); ++i)
 647         {
 648             std::fprintf(out, "#   %s\n", sc.sel[i]->selectionText());
 649         }
 650         std::fprintf(out, "#\n");
 651     }
 652 }
 653
 654 // static
 655 HelpTopicPointer
 656 SelectionCollection::createDefaultHelpTopic()
 657 {
 658     return createSelectionHelpTopic();
 659 }
 660
 661 } // namespace gmx