src/gromacs/selection/selection.h

   1 /*
   2  * This file is part of the GROMACS molecular simulation package.
   3  *
   4  * Copyright (c) 2009,2010,2011,2012,2013 by the GROMACS development team.
   5  * Copyright (c) 2014,2015,2016,2017,2018 by the GROMACS development team.
   6  * Copyright (c) 2019,2020,2021, by the GROMACS development team, led by
   7  * Mark Abraham, David van der Spoel, Berk Hess, and Erik Lindahl,
   8  * and including many others, as listed in the AUTHORS file in the
   9  * top-level source directory and at http://www.gromacs.org.
  10  *
  11  * GROMACS is free software; you can redistribute it and/or
  12  * modify it under the terms of the GNU Lesser General Public License
  13  * as published by the Free Software Foundation; either version 2.1
  14  * of the License, or (at your option) any later version.
  15  *
  16  * GROMACS is distributed in the hope that it will be useful,
  17  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  18  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  19  * Lesser General Public License for more details.
  20  *
  21  * You should have received a copy of the GNU Lesser General Public
  22  * License along with GROMACS; if not, see
  23  * http://www.gnu.org/licenses, or write to the Free Software Foundation,
  24  * Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA.
  25  *
  26  * If you want to redistribute modifications to GROMACS, please
  27  * consider that scientific software is very special. Version
  28  * control is crucial - bugs must be traceable. We will be happy to
  29  * consider code for inclusion in the official distribution, but
  30  * derived work must not be called official GROMACS. Details are found
  31  * in the README & COPYING files - if they are missing, get the
  32  * official version at http://www.gromacs.org.
  33  *
  34  * To help us fund GROMACS development, we humbly ask that you cite
  35  * the research papers on the package. Check out http://www.gromacs.org.
  36  */
  37 /*! \file
  38  * \brief
  39  * Declares gmx::Selection and supporting classes.
  40  *
  41  * \author Teemu Murtola <teemu.murtola@gmail.com>
  42  * \inpublicapi
  43  * \ingroup module_selection
  44  */
  45 #ifndef GMX_SELECTION_SELECTION_H
  46 #define GMX_SELECTION_SELECTION_H
  47
  48 #include <memory>
  49 #include <string>
  50 #include <vector>
  51
  52 #include "gromacs/utility/arrayref.h"
  53 #include "gromacs/utility/classhelpers.h"
  54 #include "gromacs/utility/gmxassert.h"
  55
  56 #include "position.h"
  57 #include "selectionenums.h"
  58
  59 struct gmx_mtop_t;
  60
  61 namespace gmx
  62 {
  63
  64 class SelectionOptionStorage;
  65 class SelectionTreeElement;
  66
  67 class AnalysisNeighborhoodPositions;
  68 class Selection;
  69 class SelectionPosition;
  70
  71 //! Container of selections used in public selection interfaces.
  72 typedef std::vector<Selection> SelectionList;
  73
  74 namespace internal
  75 {
  76
  77 /*! \internal
  78  * \brief
  79  * Internal data for a single selection.
  80  *
  81  * This class is internal to the selection module, but resides in a public
  82  * header because of efficiency reasons: it allows frequently used access
  83  * methods in \ref Selection to be inlined.
  84  *
  85  * Methods in this class do not throw unless otherwise specified.
  86  *
  87  * \ingroup module_selection
  88  */
  89 class SelectionData
  90 {
  91 public:
  92     /*! \brief
  93      * Creates a new selection object.
  94      *
  95      * \param[in] elem   Root of the evaluation tree for this selection.
  96      * \param[in] selstr String that was parsed to produce this selection.
  97      * \throws    std::bad_alloc if out of memory.
  98      */
  99     SelectionData(SelectionTreeElement* elem, const char* selstr);
 100     ~SelectionData();
 101
 102     //! Returns the name for this selection.
 103     const char* name() const { return name_.c_str(); }
 104     //! Returns the string that was parsed to produce this selection.
 105     const char* selectionText() const { return selectionText_.c_str(); }
 106     //! Returns true if the size of the selection (posCount()) is dynamic.
 107     bool isDynamic() const { return bDynamic_; }
 108     //! Returns the type of positions in the selection.
 109     e_index_t type() const { return rawPositions_.m.type; }
 110     //! Returns true if the selection only contains positions with a single atom each.
 111     bool hasOnlyAtoms() const { return type() == INDEX_ATOM; }
 112     //! Returns `true` if the atom indices in the selection are in ascending order.
 113     bool hasSortedAtomIndices() const;
 114
 115     //! Number of positions in the selection.
 116     int posCount() const { return rawPositions_.count(); }
 117     //! Returns the root of the evaluation tree for this selection.
 118     SelectionTreeElement& rootElement() { return rootElement_; }
 119
 120     //! Returns whether the covered fraction can change between frames.
 121     bool isCoveredFractionDynamic() const { return bDynamicCoveredFraction_; }
 122
 123     //! Returns true if the given flag is set.
 124     bool hasFlag(SelectionFlag flag) const { return flags_.test(flag); }
 125     //! Sets the flags for this selection.
 126     void setFlags(SelectionFlags flags) { flags_ = flags; }
 127
 128     //! \copydoc Selection::initCoveredFraction()
 129     bool initCoveredFraction(e_coverfrac_t type);
 130
 131     /*! \brief
 132      * Updates the name of the selection if missing.
 133      *
 134      * \throws    std::bad_alloc if out of memory.
 135      *
 136      * If selections get their value from a group reference that cannot be
 137      * resolved during parsing, the name is final only after group
 138      * references have been resolved.
 139      *
 140      * This function is called by SelectionCollection::setIndexGroups().
 141      */
 142     void refreshName();
 143     /*! \brief
 144      * Computes total masses and charges for all selection positions.
 145      *
 146      * \param[in] top   Topology information.
 147      * \throws    std::bad_alloc if out of memory.
 148      *
 149      * For dynamic selections, the values need to be updated after each
 150      * evaluation with refreshMassesAndCharges().
 151      * This is done by SelectionEvaluator.
 152      *
 153      * This function is called by SelectionCompiler.
 154      *
 155      * Strong exception safety guarantee.
 156      */
 157     void initializeMassesAndCharges(const gmx_mtop_t* top);
 158     /*! \brief
 159      * Updates masses and charges after dynamic selection has been
 160      * evaluated.
 161      *
 162      * \param[in] top   Topology information.
 163      *
 164      * Called by SelectionEvaluator.
 165      */
 166     void refreshMassesAndCharges(const gmx_mtop_t* top);
 167     /*! \brief
 168      * Updates the covered fraction after a selection has been evaluated.
 169      *
 170      * Called by SelectionEvaluator.
 171      */
 172     void updateCoveredFractionForFrame();
 173     /*! \brief
 174      * Computes average covered fraction after all frames have been evaluated.
 175      *
 176      * \param[in] nframes  Number of frames that have been evaluated.
 177      *
 178      * \p nframes should be equal to the number of calls to
 179      * updateCoveredFractionForFrame().
 180      * Called by SelectionEvaluator::evaluateFinal().
 181      */
 182     void computeAverageCoveredFraction(int nframes);
 183     /*! \brief
 184      * Restores position information to state it was in after compilation.
 185      *
 186      * \param[in] top   Topology information.
 187      *
 188      * Depends on SelectionCompiler storing the original atoms in the
 189      * \a rootElement_ object.
 190      * Called by SelectionEvaluator::evaluateFinal().
 191      */
 192     void restoreOriginalPositions(const gmx_mtop_t* top);
 193
 194 private:
 195     //! Name of the selection.
 196     std::string name_;
 197     //! The actual selection string.
 198     std::string selectionText_;
 199     //! Low-level representation of selected positions.
 200     gmx_ana_pos_t rawPositions_;
 201     //! Total masses for the current positions.
 202     std::vector<real> posMass_;
 203     //! Total charges for the current positions.
 204     std::vector<real> posCharge_;
 205     SelectionFlags    flags_;
 206     //! Root of the selection evaluation tree.
 207     SelectionTreeElement& rootElement_;
 208     //! Type of the covered fraction.
 209     e_coverfrac_t coveredFractionType_;
 210     //! Covered fraction of the selection for the current frame.
 211     real coveredFraction_;
 212     //! The average covered fraction (over the trajectory).
 213     real averageCoveredFraction_;
 214     //! true if the value can change as a function of time.
 215     bool bDynamic_;
 216     //! true if the covered fraction depends on the frame.
 217     bool bDynamicCoveredFraction_;
 218
 219     /*! \brief
 220      * Needed to wrap access to information.
 221      */
 222     friend class gmx::Selection;
 223     /*! \brief
 224      * Needed for proper access to position information.
 225      */
 226     friend class gmx::SelectionPosition;
 227
 228     GMX_DISALLOW_COPY_AND_ASSIGN(SelectionData);
 229 };
 230
 231 } // namespace internal
 232
 233 /*! \brief
 234  * Provides access to a single selection.
 235  *
 236  * This class provides a public interface for accessing selection information.
 237  * General information about the selection can be accessed with methods name(),
 238  * selectionText(), isDynamic(), and type().  The first three can be accessed
 239  * any time after the selection has been parsed, and type() can be accessed
 240  * after the selection has been compiled.
 241  *
 242  * There are a few methods that can be used to change the behavior of the
 243  * selection.  setEvaluateVelocities() and setEvaluateForces() can be called
 244  * before the selection is compiled to request evaluation of velocities and/or
 245  * forces in addition to coordinates.
 246  *
 247  * Each selection is made of a set of positions.  Each position has associated
 248  * coordinates, and possibly velocities and forces if they have been requested
 249  * and are available.  It also has a set of atoms associated with it; typically
 250  * the coordinates are the center-of-mass or center-of-geometry coordinates for
 251  * that set of atoms.  To access the number of positions in the selection, use
 252  * posCount().  To access individual positions, use position().
 253  * See SelectionPosition for details of how to use individual positions.
 254  * setOriginalId() can be used to adjust the return value of
 255  * SelectionPosition::mappedId(); see that method for details.
 256  *
 257  * It is also possible to access the list of atoms that make up all the
 258  * positions directly: atomCount() returns the total number of atoms in the
 259  * selection and atomIndices() an array of their indices.
 260  * Similarly, it is possible to access the coordinates and other properties
 261  * of the positions as continuous arrays through coordinates(), velocities(),
 262  * forces(), masses(), charges(), refIds(), and mappedIds().
 263  *
 264  * Both positions and atoms can be accessed after the selection has been
 265  * compiled.  For dynamic selections, the return values of these methods change
 266  * after each evaluation to reflect the situation for the current frame.
 267  * Before any frame has been evaluated, these methods return the maximal set
 268  * to which the selection can evaluate.
 269  *
 270  * There are two possible modes for how positions for dynamic selections are
 271  * handled.  In the default mode, posCount() can change, and for each frame,
 272  * only the positions that are selected in that frame can be accessed.  In a
 273  * masked mode, posCount() remains constant, i.e., the positions are always
 274  * evaluated for the maximal set, and SelectionPosition::selected() is used to
 275  * determine whether a position is selected for a frame.  The masked mode can
 276  * be requested with SelectionOption::dynamicMask().
 277  *
 278  * The class also provides methods for printing out information: printInfo()
 279  * and printDebugInfo().  These are mainly for internal use by Gromacs.
 280  *
 281  * This class works like a pointer type: copying and assignment is lightweight,
 282  * and all copies work interchangeably, accessing the same internal data.
 283  *
 284  * Methods in this class do not throw.
 285  *
 286  * \see SelectionPosition
 287  *
 288  * \inpublicapi
 289  * \ingroup module_selection
 290  */
 291 class Selection
 292 {
 293 public:
 294     /*! \brief
 295      * Creates a selection wrapper that has no associated selection.
 296      *
 297      * Any attempt to call methods in the object before a selection is
 298      * assigned results in undefined behavior.
 299      * isValid() returns `false` for the selection until it is initialized.
 300      */
 301     Selection() : sel_(nullptr) {}
 302     /*! \brief
 303      * Creates a new selection object.
 304      *
 305      * \param  sel  Selection data to wrap.
 306      *
 307      * Only for internal use by the selection module.
 308      */
 309     explicit Selection(internal::SelectionData* sel) : sel_(sel) {}
 310
 311     //! Returns whether the selection object is initialized.
 312     bool isValid() const { return sel_ != nullptr; }
 313
 314     //! Returns whether two selection objects wrap the same selection.
 315     bool operator==(const Selection& other) const { return sel_ == other.sel_; }
 316     //! Returns whether two selection objects wrap different selections.
 317     bool operator!=(const Selection& other) const { return !operator==(other); }
 318
 319     //! Returns the name of the selection.
 320     const char* name() const { return data().name(); }
 321     //! Returns the string that was parsed to produce this selection.
 322     const char* selectionText() const { return data().selectionText(); }
 323     //! Returns true if the size of the selection (posCount()) is dynamic.
 324     bool isDynamic() const { return data().isDynamic(); }
 325     //! Returns the type of positions in the selection.
 326     e_index_t type() const { return data().type(); }
 327     //! Returns true if the selection only contains positions with a single atom each.
 328     bool hasOnlyAtoms() const { return data().hasOnlyAtoms(); }
 329     //! Returns `true` if the atom indices in the selection are in ascending order.
 330     bool hasSortedAtomIndices() const { return data().hasSortedAtomIndices(); }
 331
 332     //! Total number of atoms in the selection.
 333     int atomCount() const { return data().rawPositions_.m.mapb.nra; }
 334     //! Returns atom indices of all atoms in the selection.
 335     ArrayRef<const int> atomIndices() const
 336     {
 337         return constArrayRefFromArray(sel_->rawPositions_.m.mapb.a, sel_->rawPositions_.m.mapb.nra);
 338     }
 339     //! Number of positions in the selection.
 340     int posCount() const { return data().posCount(); }
 341     //! Access a single position.
 342     SelectionPosition position(int i) const;
 343     //! Returns coordinates for this selection as a continuous array.
 344     ArrayRef<const rvec> coordinates() const
 345     {
 346         return constArrayRefFromArray(data().rawPositions_.x, posCount());
 347     }
 348     //! Returns whether velocities are available for this selection.
 349     bool hasVelocities() const { return data().rawPositions_.v != nullptr; }
 350     /*! \brief
 351      * Returns velocities for this selection as a continuous array.
 352      *
 353      * Must not be called if hasVelocities() returns false.
 354      */
 355     ArrayRef<const rvec> velocities() const
 356     {
 357         GMX_ASSERT(hasVelocities(), "Velocities accessed, but unavailable");
 358         return constArrayRefFromArray(data().rawPositions_.v, posCount());
 359     }
 360     //! Returns whether forces are available for this selection.
 361     bool hasForces() const { return sel_->rawPositions_.f != nullptr; }
 362     /*! \brief
 363      * Returns forces for this selection as a continuous array.
 364      *
 365      * Must not be called if hasForces() returns false.
 366      */
 367     ArrayRef<const rvec> forces() const
 368     {
 369         GMX_ASSERT(hasForces(), "Forces accessed, but unavailable");
 370         return constArrayRefFromArray(data().rawPositions_.f, posCount());
 371     }
 372     //! Returns masses for this selection as a continuous array.
 373     ArrayRef<const real> masses() const
 374     {
 375         // posMass_ may have more entries than posCount() in the case of
 376         // dynamic selections that don't have a topology
 377         // (and thus the masses and charges are fixed).
 378         GMX_ASSERT(data().posMass_.size() >= static_cast<size_t>(posCount()),
 379                    "Internal inconsistency");
 380         return makeArrayRef(data().posMass_).subArray(0, posCount());
 381     }
 382     //! Returns charges for this selection as a continuous array.
 383     ArrayRef<const real> charges() const
 384     {
 385         // posCharge_ may have more entries than posCount() in the case of
 386         // dynamic selections that don't have a topology
 387         // (and thus the masses and charges are fixed).
 388         GMX_ASSERT(data().posCharge_.size() >= static_cast<size_t>(posCount()),
 389                    "Internal inconsistency");
 390         return makeArrayRef(data().posCharge_).subArray(0, posCount());
 391     }
 392     /*! \brief
 393      * Returns reference IDs for this selection as a continuous array.
 394      *
 395      * \see SelectionPosition::refId()
 396      */
 397     ArrayRef<const int> refIds() const
 398     {
 399         return constArrayRefFromArray(data().rawPositions_.m.refid, posCount());
 400     }
 401     /*! \brief
 402      * Returns mapped IDs for this selection as a continuous array.
 403      *
 404      * \see SelectionPosition::mappedId()
 405      */
 406     ArrayRef<const int> mappedIds() const
 407     {
 408         return constArrayRefFromArray(data().rawPositions_.m.mapid, posCount());
 409     }
 410
 411     //! Returns whether the covered fraction can change between frames.
 412     bool isCoveredFractionDynamic() const { return data().isCoveredFractionDynamic(); }
 413     //! Returns the covered fraction for the current frame.
 414     real coveredFraction() const { return data().coveredFraction_; }
 415
 416     /*! \brief
 417      * Allows passing a selection directly to neighborhood searching.
 418      *
 419      * When initialized this way, AnalysisNeighborhoodPair objects return
 420      * indices that can be used to index the selection positions with
 421      * position().
 422      *
 423      * Works exactly like if AnalysisNeighborhoodPositions had a
 424      * constructor taking a Selection object as a parameter.
 425      * See AnalysisNeighborhoodPositions for rationale and additional
 426      * discussion.
 427      */
 428     operator AnalysisNeighborhoodPositions() const;
 429
 430     /*! \brief
 431      * Initializes information about covered fractions.
 432      *
 433      * \param[in] type Type of covered fraction required.
 434      * \returns   true if the covered fraction can be calculated for the
 435      *      selection.
 436      */
 437     bool initCoveredFraction(e_coverfrac_t type) { return data().initCoveredFraction(type); }
 438     /*! \brief
 439      * Sets whether this selection evaluates velocities for positions.
 440      *
 441      * \param[in] bEnabled  If true, velocities are evaluated.
 442      *
 443      * If you request the evaluation, but then evaluate the selection for
 444      * a frame that does not contain velocity information, results are
 445      * undefined.
 446      *
 447      * \todo
 448      * Implement it such that in the above case, hasVelocities() will
 449      * return false for such frames.
 450      *
 451      * Does not throw.
 452      */
 453     void setEvaluateVelocities(bool bEnabled)
 454     {
 455         data().flags_.set(efSelection_EvaluateVelocities, bEnabled);
 456     }
 457     /*! \brief
 458      * Sets whether this selection evaluates forces for positions.
 459      *
 460      * \param[in] bEnabled  If true, forces are evaluated.
 461      *
 462      * If you request the evaluation, but then evaluate the selection for
 463      * a frame that does not contain force information, results are
 464      * undefined.
 465      *
 466      * Does not throw.
 467      */
 468     void setEvaluateForces(bool bEnabled)
 469     {
 470         data().flags_.set(efSelection_EvaluateForces, bEnabled);
 471     }
 472
 473     /*! \brief
 474      * Sets the ID for the \p i'th position for use with
 475      * SelectionPosition::mappedId().
 476      *
 477      * \param[in] i  Zero-based index
 478      * \param[in] id Identifier to set.
 479      *
 480      * This method is not part of SelectionPosition because that interface
 481      * only provides access to const data by design.
 482      *
 483      * This method can only be called after compilation, before the
 484      * selection has been evaluated for any frame.
 485      *
 486      * \see SelectionPosition::mappedId()
 487      */
 488     void setOriginalId(int i, int id);
 489     /*! \brief
 490      * Inits the IDs for use with SelectionPosition::mappedId() for
 491      * grouping.
 492      *
 493      * \param[in] top   Topology information
 494      *     (can be NULL if not required for \p type).
 495      * \param[in] type  Type of groups to generate.
 496      * \returns   Number of groups that were present in the selection.
 497      * \throws    InconsistentInputError if the selection positions cannot
 498      *     be assigned to groups of the given type.
 499      *
 500      * If `type == INDEX_ATOM`, the IDs are initialized to 0, 1, 2, ...,
 501      * and the return value is the number of positions.
 502      * If `type == INDEX_ALL`, all the IDs are initialized to 0, and the
 503      * return value is one.
 504      * If `type == INDEX_RES` or `type == INDEX_MOL`, the first position
 505      * will get ID 0, and all following positions that belong to the same
 506      * residue/molecule will get the same ID.  The first position that
 507      * belongs to a different residue/molecule will get ID 1, and so on.
 508      * If some position contains atoms from multiple residues/molecules,
 509      * i.e., the mapping is ambiguous, an exception is thrown.
 510      * The return value is the number of residues/molecules that are
 511      * present in the selection positions.
 512      *
 513      * This method is useful if the calling code needs to group the
 514      * selection, e.g., for computing aggregate properties for each residue
 515      * or molecule.  It can then use this method to initialize the
 516      * appropriate grouping, use the return value to allocate a
 517      * sufficiently sized buffer to store the aggregated values, and then
 518      * use SelectionPosition::mappedId() to identify the location where to
 519      * aggregate to.
 520      *
 521      * \see setOriginalId()
 522      * \see SelectionPosition::mappedId()
 523      */
 524     int initOriginalIdsToGroup(const gmx_mtop_t* top, e_index_t type);
 525
 526     /*! \brief
 527      * Prints out one-line description of the selection.
 528      *
 529      * \param[in] fp      Where to print the information.
 530      *
 531      * The output contains the name of the selection, the number of atoms
 532      * and the number of positions, and indication of whether the selection
 533      * is dynamic.
 534      */
 535     void printInfo(FILE* fp) const;
 536     /*! \brief
 537      * Prints out extended information about the selection for debugging.
 538      *
 539      * \param[in] fp      Where to print the information.
 540      * \param[in] nmaxind Maximum number of values to print in lists
 541      *      (-1 = print all).
 542      */
 543     void printDebugInfo(FILE* fp, int nmaxind) const;
 544
 545 private:
 546     internal::SelectionData& data()
 547     {
 548         GMX_ASSERT(sel_ != nullptr, "Attempted to access uninitialized selection");
 549         return *sel_;
 550     }
 551     const internal::SelectionData& data() const
 552     {
 553         GMX_ASSERT(sel_ != nullptr, "Attempted to access uninitialized selection");
 554         return *sel_;
 555     }
 556
 557     /*! \brief
 558      * Pointer to internal data for the selection.
 559      *
 560      * The memory for this object is managed by a SelectionCollection
 561      * object, and the \ref Selection class simply provides a public
 562      * interface for accessing the data.
 563      */
 564     internal::SelectionData* sel_;
 565
 566     /*! \brief
 567      * Needed to access the data to adjust flags.
 568      */
 569     friend class SelectionOptionStorage;
 570 };
 571
 572 /*! \brief
 573  * Provides access to information about a single selected position.
 574  *
 575  * Each position has associated coordinates, and possibly velocities and forces
 576  * if they have been requested and are available.  It also has a set of atoms
 577  * associated with it; typically the coordinates are the center-of-mass or
 578  * center-of-geometry coordinates for that set of atoms.  It is possible that
 579  * there are not atoms associated if the selection has been provided as a fixed
 580  * position.
 581  *
 582  * After the selection has been compiled, but not yet evaluated, the contents
 583  * of the coordinate, velocity and force vectors are undefined.
 584  *
 585  * Default copy constructor and assignment operators are used, and work as
 586  * intended: the copy references the same position and works identically.
 587  *
 588  * Methods in this class do not throw.
 589  *
 590  * \see Selection
 591  *
 592  * \inpublicapi
 593  * \ingroup module_selection
 594  */
 595 class SelectionPosition
 596 {
 597 public:
 598     /*! \brief
 599      * Constructs a wrapper object for given selection position.
 600      *
 601      * \param[in] sel    Selection from which the position is wrapped.
 602      * \param[in] index  Zero-based index of the position to wrap.
 603      *
 604      * Asserts if \p index is out of range.
 605      *
 606      * Only for internal use of the library.  To obtain a SelectionPosition
 607      * object in other code, use Selection::position().
 608      */
 609     SelectionPosition(const internal::SelectionData& sel, int index) : sel_(&sel), i_(index)
 610     {
 611         GMX_ASSERT(index >= 0 && index < sel.posCount(), "Invalid selection position index");
 612     }
 613
 614     /*! \brief
 615      * Returns type of this position.
 616      *
 617      * Currently always returns the same as Selection::type().
 618      */
 619     e_index_t type() const { return sel_->type(); }
 620     //! Returns coordinates for this position.
 621     const rvec& x() const { return sel_->rawPositions_.x[i_]; }
 622     /*! \brief
 623      * Returns velocity for this position.
 624      *
 625      * Must not be called if Selection::hasVelocities() returns false.
 626      */
 627     const rvec& v() const
 628     {
 629         GMX_ASSERT(sel_->rawPositions_.v != nullptr, "Velocities accessed, but unavailable");
 630         return sel_->rawPositions_.v[i_];
 631     }
 632     /*! \brief
 633      * Returns force for this position.
 634      *
 635      * Must not be called if Selection::hasForces() returns false.
 636      */
 637     const rvec& f() const
 638     {
 639         GMX_ASSERT(sel_->rawPositions_.f != nullptr, "Forces accessed, but unavailable");
 640         return sel_->rawPositions_.f[i_];
 641     }
 642     /*! \brief
 643      * Returns total mass for this position.
 644      *
 645      * Returns the total mass of atoms that make up this position.
 646      * If there are no atoms associated or masses are not available,
 647      * returns unity.
 648      */
 649     real mass() const { return sel_->posMass_[i_]; }
 650     /*! \brief
 651      * Returns total charge for this position.
 652      *
 653      * Returns the sum of charges of atoms that make up this position.
 654      * If there are no atoms associated or charges are not available,
 655      * returns zero.
 656      */
 657     real charge() const { return sel_->posCharge_[i_]; }
 658     //! Returns the number of atoms that make up this position.
 659     int atomCount() const
 660     {
 661         return sel_->rawPositions_.m.mapb.index[i_ + 1] - sel_->rawPositions_.m.mapb.index[i_];
 662     }
 663     //! Return atom indices that make up this position.
 664     ArrayRef<const int> atomIndices() const
 665     {
 666         const int* atoms = sel_->rawPositions_.m.mapb.a;
 667         if (atoms == nullptr)
 668         {
 669             return ArrayRef<const int>();
 670         }
 671         const int first = sel_->rawPositions_.m.mapb.index[i_];
 672         return constArrayRefFromArray(&atoms[first], atomCount());
 673     }
 674     /*! \brief
 675      * Returns whether this position is selected in the current frame.
 676      *
 677      * The return value is equivalent to \c refid() == -1.  Returns always
 678      * true if SelectionOption::dynamicMask() has not been set.
 679      *
 680      * \see refId()
 681      */
 682     bool selected() const { return refId() >= 0; }
 683     /*! \brief
 684      * Returns reference ID for this position.
 685      *
 686      * For dynamic selections, this provides means to associate positions
 687      * across frames.  After compilation, these IDs are consequently
 688      * numbered starting from zero.  For each frame, the ID then reflects
 689      * the location of the position in the original array of positions.
 690      * If SelectionOption::dynamicMask() has been set for the parent
 691      * selection, the IDs for positions not present in the current
 692      * selection are set to -1, otherwise they are removed completely.
 693      *
 694      * Example:
 695      * If a dynamic selection consists of at most three positions, after
 696      * compilation refId() will return 0, 1, 2 for them, respectively.
 697      * If for a particular frame, only the first and the third are present,
 698      * refId() will return 0, 2.
 699      * If SelectionOption::dynamicMask() has been set, all three positions
 700      * can be accessed also for that frame and refId() will return 0, -1,
 701      * 2.
 702      */
 703     int refId() const { return sel_->rawPositions_.m.refid[i_]; }
 704     /*! \brief
 705      * Returns mapped ID for this position.
 706      *
 707      * Returns ID of the position that corresponds to that set with
 708      * Selection::setOriginalId().
 709      *
 710      * If for an array \c id, \c setOriginalId(i, id[i]) has been called
 711      * for each \c i, then it always holds that
 712      * \c mappedId()==id[refId()].
 713      *
 714      * Selection::setOriginalId() has not been called, the default values
 715      * are dependent on type():
 716      *  - ::INDEX_ATOM: atom indices
 717      *  - ::INDEX_RES:  residue indices
 718      *  - ::INDEX_MOL:  molecule indices
 719      *  .
 720      * All the default values are zero-based.
 721      */
 722     int mappedId() const { return sel_->rawPositions_.m.mapid[i_]; }
 723
 724     /*! \brief
 725      * Allows passing a selection position directly to neighborhood searching.
 726      *
 727      * When initialized this way, AnalysisNeighborhoodPair objects return
 728      * the index that can be used to access this position using
 729      * Selection::position().
 730      *
 731      * Works exactly like if AnalysisNeighborhoodPositions had a
 732      * constructor taking a SelectionPosition object as a parameter.
 733      * See AnalysisNeighborhoodPositions for rationale and additional
 734      * discussion.
 735      */
 736     operator AnalysisNeighborhoodPositions() const;
 737
 738 private:
 739     const internal::SelectionData* sel_;
 740     int                            i_;
 741 };
 742
 743
 744 inline SelectionPosition Selection::position(int i) const
 745 {
 746     return SelectionPosition(data(), i);
 747 }
 748
 749 } // namespace gmx
 750
 751 #endif