LOOS  v2.3.2
AtomicGroup.hpp
1 /*
2  This file is part of LOOS.
3 
4  LOOS (Lightweight Object-Oriented Structure library)
5  Copyright (c) 2008, Tod D. Romo, Alan Grossfield
6  Department of Biochemistry and Biophysics
7  School of Medicine & Dentistry, University of Rochester
8 
9  This package (LOOS) is free software: you can redistribute it and/or modify
10  it under the terms of the GNU General Public License as published by
11  the Free Software Foundation under version 3 of the License.
12 
13  This package is distributed in the hope that it will be useful,
14  but WITHOUT ANY WARRANTY; without even the implied warranty of
15  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
16  GNU General Public License for more details.
17 
18  You should have received a copy of the GNU General Public License
19  along with this program. If not, see <http://www.gnu.org/licenses/>.
20 */
21 
22 
23 
24 
25 
26 #if !defined(LOOS_ATOMICGROUP_HPP)
27 #define LOOS_ATOMICGROUP_HPP
28 
29 #include <iostream>
30 #include <string>
31 #include <stdexcept>
32 #include <vector>
33 #include <map>
34 #include <algorithm>
35 
36 #include <boost/unordered_set.hpp>
37 
38 
39 #include <loos_defs.hpp>
40 
41 
42 #include <Atom.hpp>
43 #include <XForm.hpp>
44 #include <PeriodicBox.hpp>
45 #include <utils.hpp>
46 #include <Matrix.hpp>
47 
48 #include <exceptions.hpp>
49 
50 
51 namespace loos {
52 
53 
55 
56  struct AtomSelector {
60  virtual bool operator()(const pAtom& atom) const =0;
61  virtual ~AtomSelector() { }
62  };
63 
64 
65  class AtomicGroup;
66  typedef boost::shared_ptr<AtomicGroup> pAtomicGroup;
67 
68 
70 
87  class AtomicGroup {
88  public:
89  typedef std::vector<pAtom>::iterator iterator;
90  typedef std::vector<pAtom>::const_iterator const_iterator;
91  typedef pAtom value_type;
92 
93  // Threshold for catching effectively zero singular values in
94  // the superposition code...
95  static const double superposition_zero_singular_value;
96 
97  public:
98  AtomicGroup() : _sorted(false) { }
99 
101 
104  AtomicGroup(const int n) : _sorted(true) {
105  assert(n >= 1 && "Invalid size in AtomicGroup(n)");
106  for (int i=1; i<=n; i++) {
107  pAtom pa(new Atom);
108  pa->id(i);
109  atoms.push_back(pa);
110  }
111  }
112 
115  _sorted(g._sorted),
116  atoms(g.atoms),
117  box(g.box)
118  { }
119 
120 
121  virtual ~AtomicGroup() { }
122 
124 
130  AtomicGroup copy(void) const;
131 
133 
138  virtual AtomicGroup* clone(void) const;
139 
140 
141  uint length(void) const { return(atoms.size()); }
142  uint size(void) const { return(atoms.size()); }
143  bool empty(void) const { return(atoms.empty()); }
144 
146  pAtom getAtom(const int i) const;
147 
148 #if !defined(SWIG)
149  pAtom& operator[](const int i);
151  const pAtom& operator[](const int i) const;
152 #endif
153 
155  AtomicGroup& append(pAtom pa) { atoms.push_back(pa); _sorted = false; return(*this); }
157  AtomicGroup& append(std::vector<pAtom> pas);
159  AtomicGroup& append(const AtomicGroup& grp);
160 
162  AtomicGroup& remove(pAtom pa) { deleteAtom(pa); return(*this); }
164  AtomicGroup& remove(std::vector<pAtom> pas);
166  AtomicGroup& remove(const AtomicGroup& grp);
167 
168  // Concatenation of groups and/or atoms
169  AtomicGroup& operator+=(const AtomicGroup& rhs);
170  AtomicGroup& operator+=(const pAtom& rhs);
171  AtomicGroup operator+(const AtomicGroup& rhs);
172  AtomicGroup operator+(const pAtom& rhs);
173 
175 
179  bool operator==(AtomicGroup& rhs);
180 
182  bool operator!=(AtomicGroup& rhs) {
183  return(!(operator==(rhs)));
184  }
185 
186 
187 #if !defined(SWIG)
188 
190 
193  bool operator==(const AtomicGroup& rhs) const;
194 
196  bool operator!=(const AtomicGroup& rhs) const {
197  return(!(operator==(rhs)));
198  }
199 
200 #endif // !defined(SWIG)
201 
202 
203 
205 
209  AtomicGroup subset(const int offset, const int len = 0);
210 
212  AtomicGroup excise(const int offset, const int len = 0);
213 
214 
216 
239  template<class EqualsOp> bool contains(const pAtom& p, const EqualsOp& op) {
240  const_iterator ci = std::find_if(begin(), end(), bind2nd(op, p));
241  return(ci != end());
242  }
243 
245  bool contains(const pAtom& p) { return(contains(p, AtomEquals())); }
246 
247 
249  template<class EqualsOp> bool contains(const AtomicGroup& g, const EqualsOp& op) {
250  for (const_iterator cj = g.begin(); cj != g.end(); ++cj)
251  if (std::find_if(begin(), end(), bind2nd(op, *cj)) == end())
252  return(false);
253  return(true);
254  }
255 
257  bool contains(const AtomicGroup& g) { return(contains(g, AtomEquals())); }
258 
259 
261  template<class EqualsOp> bool containsAny(const AtomicGroup& g, const EqualsOp& op)
262  {
263  for (const_iterator cj = g.begin(); cj != g.end(); ++cj)
264  if (std::find_if(begin(), end(), bind2nd(op, *cj)) != end())
265  return(true);
266  return(false);
267  }
268 
269 
271  bool containsAny(const AtomicGroup& g) { return(containsAny(g, AtomEquals())); }
272 
274 
277  template<class EqualsOp> AtomicGroup intersect(const AtomicGroup& g, const EqualsOp& op) {
278  AtomicGroup result;
279 
280  for (const_iterator cj = begin(); cj != end(); ++ cj)
281  if (std::find_if(g.begin(), g.end(), bind2nd(op, *cj)) != g.end())
282  result.addAtom(*cj);
283 
284  result.box = box;
285  return(result);
286  }
287 
289  AtomicGroup intersect(const AtomicGroup& g) { return(intersect(g, AtomEquals())); }
290 
292 
295  template<class EqualsOp> AtomicGroup merge(const AtomicGroup& g, const EqualsOp& op) {
296  AtomicGroup result = copy();
297 
298  for (const_iterator ci = g.begin(); ci != g.end(); ++ci)
299  if (std::find_if(begin(), end(), bind2nd(op, *ci)) == end())
300  result.addAtom(*ci);
301 
302  return(result);
303  }
304 
305 
307  AtomicGroup merge(const AtomicGroup& g) { return(merge(g, AtomEquals())); }
308 
309 
311  AtomicGroup select(const AtomSelector& sel) const;
312 
314 
318  std::vector<AtomicGroup> splitByUniqueSegid(void) const;
319 
321  std::vector<AtomicGroup> splitByMolecule(void) const {
322  AtomicGroup sortable = *this;
323  return(sortable.sortingSplitByMolecule());
324  }
325 
327  std::vector<AtomicGroup> splitByResidue(void) const;
328 
330  std::map<std::string, AtomicGroup> splitByName(void) const;
331 
332 
334 
341 
344 
345 
347 
354  pAtom findById(const int id) const;
355 
357  AtomicGroup groupFromID(const std::vector<int> &id_list) const;
358 
361  AtomicGroup getResidue(pAtom res);
362 
363 #if !defined(SWIG)
364  friend std::ostream& operator<<(std::ostream& os, const AtomicGroup& grp);
366 #endif
367 
368  // Some misc support routines...
369 
371  void renumber(const int start = 1, const int stride = 1);
372  int minId(void) const;
373  int maxId(void) const;
374  int minResid(void) const;
375  int maxResid(void) const;
376  int numberOfResidues(void) const;
377  int numberOfSegids(void) const;
378 
380  bool allHaveProperty(const Atom::bits& property) const;
381 
383  bool anyHaveProperty(const Atom::bits& property) const;
384 
385  // These are now deprecated in favor of the above functions...
387  bool hasBonds(void) const;
388 
390  bool hasCoords(void) const;
391 
393  void clearBonds(void);
394 
396  void pruneBonds();
397 
399  void resetAtomIndices();
400 
402  uint deduceAtomicNumberFromMass(const double tol = 0.1);
403 
405 
411  bool sorted(void) const { return(_sorted); }
412 
414  void sort(void);
415 
417  bool isPeriodic(void) const { return(box.isPeriodic()); }
418 
420  GCoord periodicBox(void) const { return(box.box()); }
421 
423  void periodicBox(const GCoord& c) { box.box(c); }
424 
426  void periodicBox(const greal x, const greal y, const greal z) {
427  box.box(GCoord(x,y,z));
428  }
429 
432 
435 
438  void reimage();
439 
441  void reimageByAtom();
442 
444  void mergeImage(pAtom &p);
446  void mergeImage();
447 
449  AtomicGroup within(const double dist, AtomicGroup& grp) const {
450  Distance2WithoutPeriodicity op;
451  return(within_private(dist, grp, op));
452  }
453 
455  AtomicGroup within(const double dist, AtomicGroup& grp, const GCoord& box) const {
456  Distance2WithPeriodicity op(box);
457  return(within_private(dist, grp, op));
458  }
459 
460 
462 
466  bool contactWith(const double dist, const AtomicGroup& grp, const uint min=1) const {
467  Distance2WithoutPeriodicity op;
468  return(contactwith_private(dist, grp, min, op));
469  }
470 
472 
476  bool contactWith(const double dist, const AtomicGroup& grp, const GCoord& box, const uint min=1) const {
477  Distance2WithPeriodicity op(box);
478  return(contactwith_private(dist, grp, min, op));
479  }
480 
481 
483 
487  // Larger distances cause problems with hydrogens...
488  void findBonds(const double dist = 1.65);
489 
490 
491 
493 
512  template<class T> T apply(T func) {
513  for (iterator i = atoms.begin(); i != atoms.end(); ++i)
514  func(*i);
515  return(func);
516  }
517 
518  // *** Helper classes...
519 
521 
537  class Iterator {
538  public:
539  explicit Iterator(const AtomicGroup& grp) : iter(grp.atoms.begin()), final(grp.atoms.end()) { }
540  pAtom operator()(void) {
541  if (iter >= final)
542  return(pAtom());
543  return(*iter++);
544  }
545  private:
546  std::vector<pAtom>::const_iterator iter, final;
547  };
548 
549  // STL-iterator access
550  // Should these reset sort status?
551  iterator begin(void) { return(atoms.begin()); }
552  iterator end(void) { return(atoms.end()); }
553 
554 #if !defined(SWIG)
555  const_iterator begin(void) const { return(atoms.begin()); }
556  const_iterator end(void) const { return(atoms.end()); }
557 #endif
558 
559  // Statistical routines...
561 
565  std::vector<GCoord> boundingBox(void) const;
566 
568 
571  GCoord centerAtOrigin(void);
572 
574  GCoord centroid(void) const;
575 
577 
581  greal radius(const bool use_atom_as_reference=false) const;
582 
584  GCoord centerOfMass(void) const;
585 
587  GCoord centerOfElectrons(void) const;
588 
590  GCoord dipoleMoment(void) const;
591  greal totalCharge(void) const;
592  greal totalMass(void) const;
593  greal radiusOfGyration(void) const;
594 
596  greal sphericalVariance(const pAtom) const;
597  greal sphericalVariance(const GCoord) const;
598 
599 
601 
604  greal rmsd(const AtomicGroup&);
605 
606  // Geometric transformations...
607 
609 
612  std::vector<GCoord> getTransformedCoords(const XForm&) const;
613 
615  void translate(const GCoord & v);
616 
618  void rotate(const GCoord& axis, const greal angle_in_degrees);
619 
621  void applyTransform(const XForm&);
622 
623 
625 
639  void copyCoordinatesFrom(const AtomicGroup& g, const uint offset = 0, const uint length = 0);
640 
642 
649  std::vector<uint> atomOrderMapFrom(const AtomicGroup& g);
650 
652 
659  void copyMappedCoordinatesFrom(const AtomicGroup& g, const std::vector<uint>& order);
660 
662 
668  void copyMappedCoordinatesFrom(const AtomicGroup& g);
669 
671  void perturbCoords(const greal);
672 
674 
699  std::vector<GCoord> principalAxes(void) const;
700 
702 
708  std::vector<GCoord> momentsOfInertia(void) const;
709 
711 
722  GMatrix superposition(const AtomicGroup&);
723 
725 
730  GMatrix alignOnto(const AtomicGroup&);
731 
732  // Set coordinates to an array
737  void setCoords(double* seq, int m, int n);
738 
739  // Return a newly allocated array containing the current group's coordinates
745  void getCoords(double** outseq, int* m, int* n);
746 
747  std::vector<double> coordsAsVector() const;
748 
749 
750  private:
751 
752 
753  // These are functors for calculating distance between two coords
754  // without and with periodicity. These can be passed to functions
755  // that need to support both ways of calculating distances, such
756  // was within_private() below...
757  struct Distance2WithoutPeriodicity {
758  double operator()(const GCoord& a, const GCoord& b) const {
759  return(a.distance2(b));
760  }
761  };
762 
763  struct Distance2WithPeriodicity {
764  Distance2WithPeriodicity(const GCoord& box) : _box(box) { }
765 
766  double operator()(const GCoord& a, const GCoord& b) const {
767  return(a.distance2(b, _box));
768  }
769 
770  GCoord _box;
771  };
772 
773 
774 
775  // Find all atoms in the current group that are within dist
776  // angstroms of any atom in the passed group. The distance
777  // calculation is determined by the passed functor so that the
778  // same code can be used for both periodic and non-periodic
779  // coordinates.
780 
781  template <typename DistanceCalc>
782  AtomicGroup within_private(const double dist, AtomicGroup& grp, const DistanceCalc& distance_functor) const {
783 
784  AtomicGroup res;
785  res.box = box;
786 
787  double dist2 = dist * dist;
788  std::vector<uint> indices;
789 
790  for (uint j=0; j<size(); j++) {
791  GCoord c = atoms[j]->coords();
792  for (uint i=0; i<grp.size(); i++) {
793  if (distance_functor(c, grp.atoms[i]->coords()) <= dist2) {
794  indices.push_back(j);
795  break;
796  }
797  }
798  }
799 
800  if (indices.size() == 0)
801  return(res);
802 
803  for (std::vector<uint>::const_iterator ci = indices.begin(); ci != indices.end(); ++ci)
804  res.addAtom(atoms[*ci]);
805 
806  return(res);
807  }
808 
809 
810  template<typename DistanceCalc>
811  bool contactwith_private(const double dist, const AtomicGroup& grp, const uint min_contacts, const DistanceCalc& distance_function) const {
812  double dist2 = dist * dist;
813  uint ncontacts = 0;
814 
815  for (uint j = 0; j<size(); ++j) {
816  GCoord c = atoms[j]->coords();
817  for (uint i = 0; i<grp.size(); ++i)
818  if (distance_function(c, grp.atoms[i]->coords()) <= dist2)
819  if (++ncontacts >= min_contacts)
820  return(true);
821  }
822  return(false);
823  }
824 
825 
826  std::vector<AtomicGroup> sortingSplitByMolecule();
827 
828  // *** Internal routines *** See the .cpp file for details...
829  void sorted(bool b) { _sorted = b; }
830 
831  pAtom findById_linearSearch(const int id) const;
832  pAtom findById_binarySearch(const int id) const;
833 
834  int rangeCheck(int) const;
835 
836  void addAtom(pAtom pa) { atoms.push_back(pa); _sorted = false; }
837  void deleteAtom(pAtom pa);
838 
839  boost::tuple<iterator, iterator> calcSubsetIterators(const int offset, const int len = 0);
840 
841  void copyCoordinatesById(AtomicGroup& g);
842 
843  // Some helper classes for using the STL
844  struct CmpById {
845  bool operator()(const pAtom& a, const pAtom& b) {
846  return(a->id() < b->id());
847  }
848  };
849 
850  struct BindId {
851  BindId(const int i) : id(i) { }
852  bool operator()(const pAtom& a) { return(a->id() == id); }
853  int id;
854  };
855 
856  typedef boost::unordered_set<int> HashInt;
857 
858  void walkBonds(AtomicGroup& mygroup, HashInt& seen, AtomicGroup& working, pAtom& moi);
859 
860 
861  double *coordsAsArray(void) const;
862  double *transformedCoordsAsArray(const XForm&) const;
863 
864  bool _sorted;
865 
866 
867  protected:
868 
869  void setGroupConnectivity();
870 
871 
872  std::vector<pAtom> atoms;
874 
875  };
876 
877  AtomicGroup operator+(const pAtom& lhs, const pAtom& rhs);
878  AtomicGroup operator+(const pAtom& lhs, const AtomicGroup& rhs);
879 
880 
881 }
882 
883 #endif
AtomicGroup subset(const int offset, const int len=0)
subset() and excise() args are patterned after perl's substr...
AtomicGroup(const AtomicGroup &g)
Copy constructor (atoms and box shared)
pAtom getAtom(const int i) const
Get the ith atom from this group.
Definition: AtomicGroup.cpp:83
bool contains(const AtomicGroup &g, const EqualsOp &op)
Determines if the passed group is a subset of the current group using the EqualsOp atom-equality poli...
bool allHaveProperty(const Atom::bits &property) const
True if all atoms in the group have the passed property(ies)
AtomicGroup select(const AtomSelector &sel) const
Return a group consisting of atoms for which sel predicate returns true...
void reimageByAtom()
Reimage atoms individually into the primary cell.
bool contains(const pAtom &p, const EqualsOp &op)
Determines if a pAtom is contained in this group using the EqualsOp atom-equality policy...
AtomicGroup intersect(const AtomicGroup &g)
Intersection of two groups.
bool operator==(AtomicGroup &rhs)
Equality test for two groups.
void findBonds(const double dist=1.65)
Distance-based search for bonds.
std::vector< AtomicGroup > splitByResidue(void) const
Returns a vector of AtomicGroups, each comprising a single residue.
AtomicGroup & append(pAtom pa)
Append the atom onto the group.
std::vector< GCoord > boundingBox(void) const
Bounding box for the group...
void periodicBox(const GCoord &c)
Set the periodic boundary conditions.
GMatrix superposition(const AtomicGroup &)
Calculates the transformation matrix for superposition of groups.
Definition: AG_linalg.cpp:188
std::map< std::string, AtomicGroup > splitByName(void) const
Returns a vector of AtomicGroups, each containing atoms with the same name.
This class manages a shared Periodicbox.
Definition: PeriodicBox.hpp:81
bool hasBonds(void) const
Does any atom in the group have bond information???
pAtom findById(const int id) const
Find a contained atom by its atomid.
AtomicGroup groupFromID(const std::vector< int > &id_list) const
Create a new group from a vector of atomids.
AtomicGroup within(const double dist, AtomicGroup &grp, const GCoord &box) const
Find atoms in grp that are within dist angstroms of atoms in the current group, considering periodici...
bool isPeriodic(void) const
Test whether or not periodic boundary conditions are set.
void renumber(const int start=1, const int stride=1)
Renumber the atomid's of the contained atoms...
AtomicGroup getResidue(pAtom res)
loos::SharedPeriodicBox sharedPeriodicBox() const
Provide access to the underlying shared periodic box...
AtomicGroup(const int n)
Creates a new AtomicGroup with n un-initialized atoms.
bool contactWith(const double dist, const AtomicGroup &grp, const uint min=1) const
Returns true if any atom of current group is within dist angstroms of grp.
Our own simple iterator for stepping over all managed atoms.
void copyCoordinatesFrom(const AtomicGroup &g, const uint offset=0, const uint length=0)
Copy coordinates from g into current group.
greal rmsd(const AtomicGroup &)
Compute the RMSD between two groups.
Basic Atom class for handling atom properties.
Definition: Atom.hpp:50
void rotate(const GCoord &axis, const greal angle_in_degrees)
Rotate group's coordinates (right-handed, about centroid)
std::vector< AtomicGroup > splitByMolecule(void) const
Returns a vector of AtomicGroups split based on bond connectivity.
AtomicGroup centrifyByResidue() const
Replace a group with the cente of masses of contained residues (see centrifyByMolecule()) ...
GCoord centerOfMass(void) const
Center of mass of the group (in group coordinates)
AtomicGroup merge(const AtomicGroup &g)
Union of two groups using the default AtomEquals atom-equality policy.
Compares two atoms based solely on name, id, resid, resname, and segid.
Definition: Atom.hpp:267
greal radius(const bool use_atom_as_reference=false) const
Maximum radius from centroid of all atoms (not gyration)
Virtual base-class for selecting atoms from a group.
Definition: AtomicGroup.hpp:56
void translate(const GCoord &v)
Translate an atomic group by vector v.
void getCoords(double **outseq, int *m, int *n)
std::vector< GCoord > momentsOfInertia(void) const
Computes the moments of inertia for a group.
Definition: AG_linalg.cpp:48
std::vector< GCoord > getTransformedCoords(const XForm &) const
Returns a vector of coordinates transformed by the passed XForm.
bits
Bits in the bitmask that flag what properties have actually been set.
Definition: Atom.hpp:65
AtomicGroup intersect(const AtomicGroup &g, const EqualsOp &op)
Computes the intersection of two groups using the EqualsOp atom-equality policy.
T apply(T func)
Apply a functor or a function to each atom in the group.
GCoord centerOfElectrons(void) const
Analogous to center of mass.
bool sorted(void) const
Is the array of atoms already sorted???
void removePeriodicBox()
Remove periodicity.
void sort(void)
Sort based on atomid.
GCoord dipoleMoment(void) const
Dipole moment, relative to group's centroid.
AtomicGroup merge(const AtomicGroup &g, const EqualsOp &op)
Union of two groups using the specified atom-equality policy.
bool contains(const pAtom &p)
Determines if a pAtom is contained in this group using the AtomEquals policy (ie the default comparis...
bool contactWith(const double dist, const AtomicGroup &grp, const GCoord &box, const uint min=1) const
Returns true if any atom of current group is within dist angstroms of grp.
GMatrix alignOnto(const AtomicGroup &)
Superimposes the current group onto the passed group.
Definition: AG_linalg.cpp:196
pAtom & operator[](const int i)
Same as getAtom(i)
Definition: AtomicGroup.cpp:90
void copyMappedCoordinatesFrom(const AtomicGroup &g, const std::vector< uint > &order)
Given a mapping of atom order, copy the coordinates into the current group.
AtomicGroup within(const double dist, AtomicGroup &grp) const
Find atoms in the current group that are within dist angstroms of any atom in grp.
virtual bool operator()(const pAtom &atom) const =0
AtomicGroup copy(void) const
Creates a deep copy of this group.
Definition: AtomicGroup.cpp:56
uint deduceAtomicNumberFromMass(const double tol=0.1)
Deduce atomic number from mass (if present), returning number of atoms assigned.
std::vector< uint > atomOrderMapFrom(const AtomicGroup &g)
Map the order of atoms in AtomicGroup g into the current group.
void clearBonds(void)
Remove any bonding information present in contained atoms.
bool operator!=(AtomicGroup &rhs)
Inequality test for two groups.
friend std::ostream & operator<<(std::ostream &os, const AtomicGroup &grp)
Output the group in pseudo-XML format...
AtomicGroup excise(const int offset, const int len=0)
excise returns the excised atoms as a group...
bool operator!=(const AtomicGroup &rhs) const
Inequality test for two groups.
virtual AtomicGroup * clone(void) const
Creates a lightweight clone of this group (for polymorphism)
Definition: AtomicGroup.cpp:51
Class for handling groups of Atoms (pAtoms, actually)
Definition: AtomicGroup.hpp:87
void pruneBonds()
Attempt to prune connectivity (only retain bonds to atoms within this AtomicGroup) ...
void applyTransform(const XForm &)
Apply the given transform to the group's coordinates...
void mergeImage()
Takes a group that's split across a periodic boundary and reimages it so it's all together...
Namespace for most things not already encapsulated within a class.
bool anyHaveProperty(const Atom::bits &property) const
True if any atom in the group have the passed property(ies)
std::vector< GCoord > principalAxes(void) const
Compute the principal axes of a group.
Definition: AG_linalg.cpp:104
bool hasCoords(void) const
Do all the atoms in the group have coordinates?
bool containsAny(const AtomicGroup &g)
Determines if a group contains any atom using the default AtomEquals policy.
void resetAtomIndices()
Reset the atom indices (used for interfacing with trajectories)
bool contains(const AtomicGroup &g)
Determines if a group is a subset of the current group using the default AtomEquals policy...
void perturbCoords(const greal)
Each atom is moved in a random direction by a vector of the passed size.
std::vector< AtomicGroup > splitByUniqueSegid(void) const
Returns a vector of AtomicGroups split from the current group based on segid.
GCoord centerAtOrigin(void)
Translates the group so that the centroid is at the origin.
void setCoords(double *seq, int m, int n)
bool containsAny(const AtomicGroup &g, const EqualsOp &op)
Determines if a group contains any atom.
AtomicGroup centrifyByMolecule() const
Replace a group with the center of masses of contained molecules.
GCoord centroid(void) const
Centroid of atoms (ignores mass, operates in group coordinates)
GCoord periodicBox(void) const
Fetch the periodic boundary conditions.
greal sphericalVariance(const pAtom) const
Spherical variance of group with respect to target atom.
void periodicBox(const greal x, const greal y, const greal z)
Set the periodic boundary conditions.