Pteros  2.0
Molecular modeling library for human beings!
selection.h
1 /*
2  *
3  * This source code is part of
4  * ******************
5  * *** Pteros ***
6  * ******************
7  * molecular modeling library
8  *
9  * Copyright (c) 2009-2015, Semen Yesylevskyy
10  *
11  * All works, which use Pteros, should cite:
12  *
13  * Semen O. Yesylevskyy, "Pteros: Fast and easy to use open-source
14  * C++ library for molecular analysis",
15  * Journal of Computational Chemistry, 2012, 33(19), 1632–1636.
16  *
17  * This program is free software; you can redistribute it and/or
18  * modify it under the terms of Artistic License:
19  * Read http://www.opensource.org/licenses/artistic-license-2.0.php
20  *
21 */
22 
23 #ifndef SELECTION_H
24 #define SELECTION_H
25 
26 #include <iostream>
27 
28 #include <string>
29 #include <memory>
30 #include <map>
31 #include <functional>
32 #include <Eigen/Core>
33 #include <Eigen/Geometry>
34 #include "pteros/core/system.h"
35 #include "pteros/core/typedefs.h"
36 
37 namespace pteros {
38 
39 // Forward declaration of friend classes
40 class System;
41 class Selection_parser;
42 class Atom_proxy;
43 
58 class Selection {
59  // System and Selection are friends because they are closely integrated.
60  friend class System;
61  friend class Selection_parser;
62  friend class Grid_searcher;
63 
64  public:
65  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
68 
70  Selection();
71 
75  Selection(const System& sys);
76 
81  Selection(const System& sys, std::string str, int fr = 0);
82 
91  Selection(const System& sys, int ind1, int ind2);
92 
95  Selection(const System& sys, const std::vector<int>& ind);
96 
99  Selection(const System& sys,
100  std::vector<int>::iterator it1,
101  std::vector<int>::iterator it2);
102 
111  Selection(const System& sys,
112  const std::function<void(const System&,int,std::vector<int>&)>& callback,
113  int fr = 0);
114 
116  Selection(const Selection& sel);
117 
119  virtual ~Selection();
120 
123 
126  bool operator==(const Selection &other) const;
127 
130  bool operator!=(const Selection &other) const {
131  return !(*this == other);
132  }
133 
155  Atom_proxy operator[](int ind) const;
156 
159  friend std::ostream& operator<<(std::ostream& os, const Selection& sel);
160 
163  friend Selection operator|(const Selection& sel1, const Selection& sel2);
164 
167  friend Selection operator&(const Selection& sel1, const Selection& sel2);
168 
172  friend Selection operator-(const Selection& sel1, const Selection& sel2);
173 
176  Selection operator~() const;
178 
179  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
182 
184  void append(const Selection& sel);
185 
187  void append(int ind);
188 
190  void remove(const Selection& sel);
191 
193  void remove(int ind);
194 
196  void invert();
197 
200  void set_system(const System& sys);
201 
206  void modify(std::string str, int fr = 0);
207 
209  void modify(int ind1, int ind2);
210 
213  void modify(const std::vector<int>& ind);
214 
217  void modify(std::vector<int>::iterator it1, std::vector<int>::iterator it2);
218 
227  void modify(const std::function<void(const System&,int,std::vector<int>&)>& callback, int fr = 0);
228 
230  void modify(const System& sys, std::string str, int fr = 0);
231 
233  void modify(const System& sys, int ind1, int ind2);
234 
236  void modify(const System& sys, const std::vector<int>& ind);
237 
239  void modify(const System& sys,
240  std::vector<int>::iterator it1,
241  std::vector<int>::iterator it2);
242 
244  void modify(const System& sys, const std::function<void(const System&,int,std::vector<int>&)>& callback, int fr = 0);
245 
251  void apply();
252 
258  void update();
259 
264  void clear();
266 
267  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
269 
299  Selection select(std::string str);
301  Selection operator()(std::string str);
302 
304  Selection select(int ind1, int ind2);
305  Selection operator()(int ind1, int ind2);
306 
308  Selection select(const std::vector<int>& ind);
309  Selection operator()(const std::vector<int>& ind);
311 
312  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
315 
317  std::vector<int>::const_iterator index_begin() const { return index.begin(); }
318 
320  std::vector<int>::const_iterator index_end() const { return index.end(); }
321 
323  std::back_insert_iterator<std::vector<int> > index_back_inserter() {
324  return std::back_inserter(index);
325  }
326 
329  class iterator;
330 
332  iterator begin();
334  iterator end();
336 
337 
338  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
341 
343  int get_frame() const {return frame;}
344 
347  void set_frame(int fr);
348 
350  System* get_system() const { return system; }
351 
354  std::string get_text() const;
355 
357  std::vector<int> get_index() const { return index; }
358 
360  std::vector<char> get_chain() const;
361 
364  void set_chain(const std::vector<char>& data);
365 
367  void set_chain(char data);
368 
370  std::vector<char> get_unique_chain() const;
371 
374  std::vector<int> get_resid() const;
375 
377  std::vector<int> get_unique_resid() const;
378 
380  std::vector<std::string> get_unique_resname() const;
381 
384  void set_resid(const std::vector<int>& data);
385 
387  void set_resid(int data);
388 
391  std::vector<int> get_resindex() const;
392 
394  std::vector<int> get_unique_resindex() const;
395 
397  std::vector<std::string> get_name() const;
398 
401  void set_name(const std::vector<std::string>& data);
402 
404  void set_name(std::string& data);
405 
407  std::vector<std::string> get_resname() const;
408 
411  void set_resname(const std::vector<std::string>& data);
412 
414  void set_resname(std::string& data);
415 
417  Eigen::MatrixXf get_xyz() const;
418 
420  void get_xyz(MatrixXf_ref res) const;
421 
423  void set_xyz(MatrixXf_const_ref coord);
424 
426  std::vector<float> get_mass() const;
427 
430  void set_mass(const std::vector<float> m);
431 
433  void set_mass(float data);
434 
436  std::vector<float> get_beta() const;
437 
440  void set_beta(std::vector<float>& data);
441 
443  void set_beta(float data);
444 
446  std::vector<float> get_occupancy() const;
447 
450  void set_occupancy(std::vector<float>& data);
451 
453  void set_occupancy(float data);
454 
456  std::vector<std::string> get_tag() const;
457 
460  void set_tag(std::vector<std::string>& data);
461 
463  void set_tag(std::string data);
465 
466 
467  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
470 
481  Eigen::Vector3f center(bool mass_weighted = false, bool periodic = false) const;
482 
484  void minmax(Vector3f_ref min, Vector3f_ref max) const;
485 
487  float sasa(float probe_r = 0.14, float* total_volume = nullptr,
488  std::vector<float>* area_per_atom = nullptr,
489  std::vector<float>* volume_per_atom = nullptr) const;
490 
492  Eigen::MatrixXf average_structure(int b=0, int e=-1) const;
493 
498  Eigen::MatrixXf atom_traj(int ind, int b=0, int e=-1) const;
499 
507  void inertia(Vector3f_ref moments, Matrix3f_ref axes, bool periodic = false) const;
508 
516  float gyration(bool periodic = false) const;
517 
521  float distance(int i, int j, bool is_periodic = true,
522  Vector3i_const_ref dims = Eigen::Vector3i::Ones()) const;
523 
527  float angle(int i, int j, int k, bool is_periodic = true,
528  Vector3i_const_ref dims = Eigen::Vector3i::Ones()) const;
529 
533  float dihedral(int i, int j, int k, int l, bool is_periodic = true,
534  Vector3i_const_ref dims = Eigen::Vector3i::Ones()) const;
535 
537 
538 
539  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
542 
544  void translate(Vector3f_const_ref v);
545 
549  void rotate(int axis, float angle);
550 
555  void rotate(int axis, float angle, Vector3f_const_ref pivot);
556 
561  void rotate(Vector3f_const_ref direction, float angle, Vector3f_const_ref pivot);
562 
564  void rotate(Matrix3f_const_ref m);
565 
567  void rotate(Vector3f_const_ref angles, Vector3f_const_ref pivot);
568 
570  void wrap(Vector3i_const_ref dims = Eigen::Vector3i::Ones());
571 
578  void unwrap(Vector3i_const_ref dims = Eigen::Vector3i::Ones());
579 
587  int unwrap_bonds(float d = 0.2, int leading_index = 0,
588  Vector3i_const_ref dims = Eigen::Vector3i::Ones());
589 
595  Eigen::Affine3f principal_transform(bool is_periodic = false) const;
596 
604  void principal_orient(bool is_periodic = false);
606 
607 
608  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
611 
613  float rmsd(int fr) const;
614 
616  float rmsd(int fr1, int fr2) const;
617 
619  friend float rmsd(const Selection& sel1, const Selection& sel2);
620 
627  friend float rmsd(const Selection& sel1, int fr1, const Selection& sel2, int fr2);
628 
630  friend void fit(Selection& sel1, const Selection& sel2);
631 
633  void fit_trajectory(int ref_frame=0, int b=0, int e=-1);
634 
636  friend Eigen::Affine3f fit_transform(const Selection& sel1, const Selection& sel2);
637 
639  Eigen::Affine3f fit_transform(int fr1, int fr2) const;
640 
642  void fit(int fr1, int fr2);
643 
645  void apply_transform(const Eigen::Affine3f& t);
647 
648 
649  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
652 
656  Energy_components non_bond_energy(float cutoff=0.25, bool periodic=true) const;
657 
661  // Default value should be given in definition of friend function!
662  friend Energy_components non_bond_energy(const Selection& sel1,
663  const Selection& sel2,
664  float cutoff,
665  int fr,
666  bool periodic);
667 
669 
670 
671  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
674 
680  // Can't be made const because of internal calls
681  void write(std::string fname, int b=-1,int e=-1);
682 
683  void write(const std::unique_ptr<Mol_file>& handler, Mol_file_content what,int b=-1,int e=-1);
685 
686 
687  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
690 
692  int size() const {return index.size();}
693 
696  bool text_based() const {
697  return sel_text!="";
698  }
699 
702  bool coord_dependent() const {
703  return (bool)parser;
704  }
705 
709  void flatten();
710 
713  std::string gromacs_ndx(std::string name);
714 
716 
717  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
720 
724  void split_by_connectivity(float d, std::vector<Selection>& res, bool periodic=true);
725 
727  void split_by_residue(std::vector<Selection>& res);
728 
730  void split_by_chain(std::vector<Selection>& chains);
731 
733  void split_by_contiguous_index(std::vector<Selection>& parts);
734 
736  void split_by_contiguous_residue(std::vector<Selection>& parts);
737 
742  template<class F>
743  void split(std::vector<Selection>& parts, F callback){
744  using namespace std;
745  parts.clear();
746  using Ret = decltype(callback(*this,0));
747  // Map
748  map<Ret,vector<int> > m;
749  for(int i=0; i<size(); ++i){
750  m[callback(*this,i)].push_back(Index(i));
751  }
752  // Create selections
753  typename map<Ret,vector<int> >::iterator it;
754  for(it=m.begin();it!=m.end();it++){
755  parts.push_back(Selection(*system));
756  parts.back().modify( it->second.begin(), it->second.end() );
757  }
758  }
759 
760 
764  void each_residue(std::vector<Selection>& sel) const;
766 
767 
768  //~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
778 
781  inline float& X(int ind){
782  return system->traj[frame].coord[index[ind]](0);
783  }
784 
785  inline const float& X(int ind) const {
786  return system->traj[frame].coord[index[ind]](0);
787  }
788 
790  inline float& X(int ind, int fr){
791  return system->traj[fr].coord[index[ind]](0);
792  }
793 
794  inline const float& X(int ind, int fr) const {
795  return system->traj[fr].coord[index[ind]](0);
796  }
797 
799  inline float& Y(int ind){
800  return system->traj[frame].coord[index[ind]](1);
801  }
802 
803  inline const float& Y(int ind) const {
804  return system->traj[frame].coord[index[ind]](1);
805  }
806 
808  inline float& Y(int ind, int fr){
809  return system->traj[fr].coord[index[ind]](1);
810  }
811 
812  inline const float& Y(int ind, int fr) const {
813  return system->traj[fr].coord[index[ind]](1);
814  }
815 
817  inline float& Z(int ind){
818  return system->traj[frame].coord[index[ind]](2);
819  }
820 
821  inline const float& Z(int ind) const {
822  return system->traj[frame].coord[index[ind]](2);
823  }
824 
826  inline float& Z(int ind, int fr){
827  return system->traj[fr].coord[index[ind]](2);
828  }
829 
830  inline const float& Z(int ind, int fr) const {
831  return system->traj[fr].coord[index[ind]](2);
832  }
833 
835  inline Eigen::Vector3f& XYZ(int ind){
836  return system->traj[frame].coord[index[ind]];
837  }
838 
839  inline const Eigen::Vector3f& XYZ(int ind) const {
840  return system->traj[frame].coord[index[ind]];
841  }
842 
845  inline Eigen::Vector3f* XYZ_ptr(int ind) const {
846  return &(system->traj[frame].coord[index[ind]]);
847  }
848 
850  inline Eigen::Vector3f& XYZ(int ind, int fr){
851  return system->traj[fr].coord[index[ind]];
852  }
853 
854  inline const Eigen::Vector3f& XYZ(int ind, int fr) const {
855  return system->traj[fr].coord[index[ind]];
856  }
857 
859  inline int& Type(int ind){
860  return system->atoms[index[ind]].type;
861  }
862 
863  inline const int& Type(int ind) const {
864  return system->atoms[index[ind]].type;
865  }
866 
868  inline std::string& Type_name(int ind){
869  return system->atoms[index[ind]].type_name;
870  }
871 
872  inline const std::string& Type_name(int ind) const {
873  return system->atoms[index[ind]].type_name;
874  }
875 
877  inline std::string& Resname(int ind){
878  return system->atoms[index[ind]].resname;
879  }
880 
881  inline const std::string& Resname(int ind) const {
882  return system->atoms[index[ind]].resname;
883  }
884 
886  inline char& Chain(int ind){
887  return system->atoms[index[ind]].chain;
888  }
889 
890  inline const char& Chain(int ind) const {
891  return system->atoms[index[ind]].chain;
892  }
893 
895  inline std::string& Name(int ind){
896  return system->atoms[index[ind]].name;
897  }
898 
899  inline const std::string& Name(int ind) const {
900  return system->atoms[index[ind]].name;
901  }
902 
904  inline float& Mass(int ind){
905  return system->atoms[index[ind]].mass;
906  }
907 
908  inline const float& Mass(int ind) const {
909  return system->atoms[index[ind]].mass;
910  }
911 
913  inline float& Charge(int ind){
914  return system->atoms[index[ind]].charge;
915  }
916 
917  inline const float& Charge(int ind) const {
918  return system->atoms[index[ind]].charge;
919  }
920 
922  inline float& Beta(int ind){
923  return system->atoms[index[ind]].beta;
924  }
925 
926  inline const float& Beta(int ind) const {
927  return system->atoms[index[ind]].beta;
928  }
929 
931  inline float& Occupancy(int ind){
932  return system->atoms[index[ind]].occupancy;
933  }
934 
935  inline const float& Occupancy(int ind) const {
936  return system->atoms[index[ind]].occupancy;
937  }
938 
940  inline int& Resid(int ind){
941  return system->atoms[index[ind]].resid;
942  }
943 
944  inline const int& Resid(int ind) const {
945  return system->atoms[index[ind]].resid;
946  }
947 
949  inline int& Index(int ind){
950  return index[ind];
951  }
952 
953  inline const int& Index(int ind) const {
954  return index[ind];
955  }
956 
958  inline std::string& Tag(int ind){
959  return system->atoms[index[ind]].tag;
960  }
961 
962  inline const std::string& Tag(int ind) const {
963  return system->atoms[index[ind]].tag;
964  }
965 
967  inline pteros::Atom& Atom_data(int ind){
968  return system->atoms[index[ind]];
969  }
970 
971  inline const pteros::Atom& Atom_data(int ind) const {
972  return system->atoms[index[ind]];
973  }
974 
976  inline int& Resindex(int ind){
977  return system->atoms[index[ind]].resindex;
978  }
979 
980  inline const int& Resindex(int ind) const {
981  return system->atoms[index[ind]].resindex;
982  }
983 
985  inline float VDW(int ind) const {
986  switch(system->atoms[index[ind]].name[0]){
987  case 'H': return 0.1;
988  case 'C': return 0.17;
989  case 'N': return 0.1625;
990  case 'O': return 0.149; //mean value used
991  case 'S': return 0.1782;
992  case 'P': return 0.1871;
993  default: return 0.17;
994  }
995  }
996 
1005  inline Periodic_box& Box() {
1006  return system->traj[frame].box;
1007  }
1008 
1009  inline const Periodic_box& Box() const {
1010  return system->traj[frame].box;
1011  }
1012 
1021  inline float& Time() {
1022  return system->traj[frame].time;
1023  }
1024 
1025  inline const float& Time() const {
1026  return system->traj[frame].time;
1027  }
1028 
1030 
1031 protected:
1032  // Row text of selection
1033  std::string sel_text;
1034  // Indexes of atoms in selection
1035  std::vector<int> index;
1036  // Pointer to target system
1037  System* system;
1038 
1039  // Stores current frame
1040  int frame;
1041 
1042  // Holds an instance of selection parser
1043  std::unique_ptr<Selection_parser> parser;
1044  void allocate_parser();
1045  void sort_and_remove_duplicates();
1046 };
1047 
1048 //==============================================================================
1049 
1053 class Atom_proxy {
1054  friend class Selection::iterator;
1055 public:
1056  Atom_proxy(){}
1057  Atom_proxy(Selection* s, int i): sel(s), ind(i) {}
1058 
1061  inline int& Resid(){ return sel->Resid(ind); }
1062  inline const int& Resid() const { return sel->Resid(ind); }
1063 
1064  inline int& Index(){ return sel->Index(ind); }
1065  inline const int& Index() const { return sel->Index(ind); }
1066 
1067  inline std::string& Name(){ return sel->Name(ind); }
1068  inline const std::string& Name() const { return sel->Name(ind); }
1069 
1070  inline char& Chain(){ return sel->Chain(ind); }
1071  inline const char& Chain() const { return sel->Chain(ind); }
1072 
1073  inline std::string& Resname(){ return sel->Resname(ind); }
1074  inline const std::string& Resname() const { return sel->Resname(ind); }
1075 
1076  inline std::string& Tag(){ return sel->Tag(ind); }
1077  inline const std::string& Tag() const { return sel->Tag(ind); }
1078 
1079  inline float& Occupancy(){ return sel->Occupancy(ind); }
1080  inline const float& Occupancy() const { return sel->Occupancy(ind); }
1081 
1082  inline float& Beta(){ return sel->Beta(ind); }
1083  inline const float& Beta() const { return sel->Beta(ind); }
1084 
1085  inline int& Resindex(){ return sel->Resindex(ind); }
1086  inline const int& Resindex() const { return sel->Resindex(ind); }
1087 
1088  inline float& Mass(){ return sel->Mass(ind); }
1089  inline const float& Mass() const { return sel->Mass(ind); }
1090 
1091  inline float& Charge(){ return sel->Charge(ind); }
1092  inline const float& Charge() const { return sel->Charge(ind); }
1093 
1094  inline int& Type(){ return sel->Type(ind); }
1095  inline const int& Type() const { return sel->Type(ind); }
1096 
1097  inline std::string& Type_name(){ return sel->Type_name(ind); }
1098  inline const std::string& Type_name() const { return sel->Type_name(ind); }
1099 
1100  inline float& X(){ return sel->X(ind); }
1101  inline const float& X() const { return sel->X(ind); }
1102 
1103  inline float& Y(){ return sel->Y(ind); }
1104  inline const float& Y() const { return sel->Y(ind); }
1105 
1106  inline float& Z(){ return sel->Z(ind); }
1107  inline const float& Z() const { return sel->Z(ind); }
1108 
1109  inline Eigen::Vector3f& XYZ(){ return sel->XYZ(ind); }
1110  inline const Eigen::Vector3f& XYZ() const { return sel->XYZ(ind); }
1111 
1112  inline float& X(int fr){ return sel->X(ind,fr); }
1113  inline const float& X(int fr) const { return sel->X(ind,fr); }
1114 
1115  inline float& Y(int fr){ return sel->Y(ind,fr); }
1116  inline const float& Y(int fr) const { return sel->Y(ind,fr); }
1117 
1118  inline float& Z(int fr){ return sel->Z(ind,fr); }
1119  inline const float& Z(int fr) const { return sel->Z(ind,fr); }
1120 
1121  inline Eigen::Vector3f& XYZ(int fr){ return sel->XYZ(ind,fr); }
1122  inline const Eigen::Vector3f& XYZ(int fr) const { return sel->XYZ(ind,fr); }
1123 
1124  inline Atom& Atom_data(){ return sel->Atom_data(ind); }
1125  inline const Atom& Atom_data() const { return sel->Atom_data(ind); }
1126 
1128 
1130  bool operator==(const Atom_proxy& other) const {
1131  return (sel==other.sel && ind==other.ind);
1132  }
1133 
1135  bool operator!=(const Atom_proxy &other) const {
1136  return !(*this == other);
1137  }
1138 
1139 protected:
1140  Selection* sel;
1141  int ind;
1142 };
1143 
1144 //==============================================================================
1145 
1147 class Selection::iterator {
1148 public:
1149  typedef Atom_proxy value_type;
1150  typedef int difference_type;
1151  typedef Atom_proxy* pointer;
1152  typedef Atom_proxy& reference;
1153  typedef std::forward_iterator_tag iterator_category;
1154 
1155  iterator(Selection* sel, int pos) { proxy.sel = sel; proxy.ind = pos; }
1156  iterator operator++() { iterator tmp = *this; proxy.ind++; return tmp; }
1157  iterator operator++(int junk) { proxy.ind++; return *this; }
1158  Atom_proxy& operator*() { return proxy; }
1159  Atom_proxy* operator->() { return &proxy; }
1160  bool operator==(const iterator& rhs) { return proxy == rhs.proxy; }
1161  bool operator!=(const iterator& rhs) { return proxy != rhs.proxy; }
1162 private:
1163  Atom_proxy proxy;
1164 };
1165 
1166 } // namespace pteros
1167 #endif /* SELECTION_H */
Eigen::Vector3f & XYZ(int ind)
Extracts X,Y and Z for current frame.
Definition: selection.h:835
float dihedral(int i, int j, int k, int l, bool is_periodic=true, Vector3i_const_ref dims=Eigen::Vector3i::Ones()) const
Get dihedral angle in degrees between three atoms (periodic in given dimensions if needed)...
Definition: selection.cpp:1862
friend Selection operator|(const Selection &sel1, const Selection &sel2)
Creates new Selection, which is the logical OR of two parent selections.
Definition: selection.cpp:504
Periodic_box & Box()
Returns periodic box of the frame pointed by selection The same as:
Definition: selection.h:1005
float rmsd(int fr) const
RMSD between current and another frame.
Definition: selection.cpp:1189
void principal_orient(bool is_periodic=false)
Orient molecule by its principal axes.
Definition: selection.cpp:1988
System * get_system() const
Get pointer to the system, which is pointed by this selection.
Definition: selection.h:350
Selection & operator=(Selection sel)
Assignment operator.
Definition: selection.cpp:450
float gyration(bool periodic=false) const
Computes radius of gyration for selection.
Definition: selection.cpp:1833
friend Selection operator-(const Selection &sel1, const Selection &sel2)
Creates new Selection, by removing all atoms of sel2 from sel1.
Definition: selection.cpp:549
friend Selection operator&(const Selection &sel1, const Selection &sel2)
Creates new Selection, which is the logical AND of two parent selections.
Definition: selection.cpp:522
Eigen::Vector3f * XYZ_ptr(int ind) const
Returns pointer to the coordinates of given atom for current frame.
Definition: selection.h:845
int & Resid(int ind)
Extracts residue number.
Definition: selection.h:940
std::vector< char > get_chain() const
Get vector of all chains in selection.
Definition: selection.cpp:641
std::vector< int >::const_iterator index_end() const
Get const iterator for the end of index.
Definition: selection.h:320
iterator begin()
Begin iterator.
Definition: selection.cpp:617
std::vector< float > get_occupancy() const
Get occupancy.
Definition: selection.cpp:911
int size() const
Get the size of selection.
Definition: selection.h:692
Pteros namespace.
Definition: options.cpp:32
std::vector< std::string > get_unique_resname() const
Get vector of unique resnames in selection.
Definition: selection.cpp:736
Atom_proxy operator[](int ind) const
Indexing operator.
Definition: selection.cpp:478
char & Chain(int ind)
Extracts chain.
Definition: selection.h:886
bool operator!=(const Selection &other) const
Inequality operator Selection are compared by their indexes.
Definition: selection.h:130
std::vector< std::string > get_tag() const
Get tags.
Definition: selection.cpp:939
void split_by_contiguous_residue(std::vector< Selection > &parts)
Split selection into contiguous ranges of resindexes.
Definition: selection.cpp:1756
float & Y(int ind)
Extracts Y for current frame.
Definition: selection.h:799
Auxilary type used to incapsulate the atom and its current coordinates Used internally in Selection::...
Definition: atom_proxy.h:1012
The system of atoms.
Definition: system.h:95
void flatten()
"Flattens" selection by removing coordinate dependence and making it not text-based.
Definition: selection.cpp:1551
iterator end()
End iterator.
Definition: selection.cpp:621
STL namespace.
bool text_based() const
Returns true if selection was created from text string and false if it was constructed &#39;by hand&#39; by a...
Definition: selection.h:696
void set_resname(const std::vector< std::string > &data)
Set resnames in selection from supplied vector.
Definition: selection.cpp:815
std::string get_text() const
Get selection text If selection is not textual returns generated index-based text "index i j k...
Definition: selection.cpp:629
float angle(int i, int j, int k, bool is_periodic=true, Vector3i_const_ref dims=Eigen::Vector3i::Ones()) const
Get angle in degrees between three atoms (periodic in given dimensions if needed).
Definition: selection.cpp:1857
std::vector< int > get_unique_resindex() const
Get vector of unique resindexes&#39;s in selection.
Definition: selection.cpp:747
std::vector< char > get_unique_chain() const
Get vector of unique chains in selection.
Definition: selection.cpp:669
pteros::Atom & Atom_data(int ind)
Extracts whole atom.
Definition: selection.h:967
void set_system(const System &sys)
Sets new system for selection.
Definition: selection.cpp:262
void split_by_residue(std::vector< Selection > &res)
Split selection by residue index.
Definition: selection.cpp:1707
Energy_components non_bond_energy(float cutoff=0.25, bool periodic=true) const
Self-energy of selection computed within given interaction cut-off.
Definition: selection.cpp:1208
Selection class.
Definition: atom_proxy.h:57
void split_by_connectivity(float d, std::vector< Selection > &res, bool periodic=true)
Split current selection into several selections according to the interatomic distances.
Definition: selection.cpp:1622
Components of the non-bond energy.
Definition: system.h:40
friend void fit(Selection &sel1, const Selection &sel2)
Fit two selection of the same size. sel1 is modified to be fit to sel2.
Definition: selection.cpp:1407
std::back_insert_iterator< std::vector< int > > index_back_inserter()
Get back_insert_iterator for index.
Definition: selection.h:323
Selection select(std::string str)
Sub-selections allow selecting atoms inside existing selection (narrowing or refining existing select...
Definition: selection.cpp:279
std::vector< int > get_index() const
Get vector of all indexes in selection.
Definition: selection.h:357
int get_frame() const
Get current frame selection is pointing to.
Definition: selection.h:343
std::string & Type_name(int ind)
Extracts typename.
Definition: selection.h:868
Eigen::MatrixXf average_structure(int b=0, int e=-1) const
Computes average structure over the range of frames.
Definition: selection.cpp:851
Eigen::Affine3f principal_transform(bool is_periodic=false) const
Get transform for orienting selection by principal axes.
Definition: selection.cpp:1953
float & Z(int ind, int fr)
Extracts Z for given frame frame fr.
Definition: selection.h:826
float & X(int ind, int fr)
Extracts X for given frame frame fr.
Definition: selection.h:790
Class which represents a single atom.
Definition: atom.h:33
void modify(std::string str, int fr=0)
Modifies selection string in existing selection.
Definition: selection.cpp:354
void set_resid(const std::vector< int > &data)
Set resid&#39;s in selection from supplied vector.
Definition: selection.cpp:689
std::string & Tag(int ind)
Extracts tag.
Definition: selection.h:958
void split_by_contiguous_index(std::vector< Selection > &parts)
Split selection into contiguous ranges of indexes.
Definition: selection.cpp:1741
void set_tag(std::vector< std::string > &data)
Set tags in selection to the values from supplied vector.
Definition: selection.cpp:949
void invert()
Inverts selection in place by selection those atoms which were not selected.
Definition: selection.cpp:238
void inertia(Vector3f_ref moments, Matrix3f_ref axes, bool periodic=false) const
Computes the central momens of inertia and principal axes of inertia.
Definition: selection.cpp:1771
std::vector< float > get_beta() const
Get beta.
Definition: selection.cpp:882
void fit_trajectory(int ref_frame=0, int b=0, int e=-1)
Fit all frames in the trajectory to reference frame.
Definition: selection.cpp:1415
Selection()
Default constructor for absolutely empty selection.
Definition: selection.cpp:78
void wrap(Vector3i_const_ref dims=Eigen::Vector3i::Ones())
Wraps whole selection to the periodic box.
Definition: selection.cpp:1867
Random-access forward iterator for Selection.
Definition: atom_proxy.h:1106
void translate(Vector3f_const_ref v)
Translate selection by given vector.
Definition: selection.cpp:1058
void append(const Selection &sel)
Append another selection to this one.
Definition: selection.cpp:193
Eigen::MatrixXf atom_traj(int ind, int b=0, int e=-1) const
Extracts X,Y,Z for given atom index for specified range of frames (gets trajectory of given atom)...
Definition: selection.cpp:1607
void write(std::string fname, int b=-1, int e=-1)
Write structure or trajectory for selection.
Definition: selection.cpp:1500
float VDW(int ind) const
Computes VDW radius. Read only.
Definition: selection.h:985
int unwrap_bonds(float d=0.2, int leading_index=0, Vector3i_const_ref dims=Eigen::Vector3i::Ones())
Unwraps selection to make it whole (without jumps over periodic box boundary).
Definition: selection.cpp:1880
float & Mass(int ind)
Extracts atom mass.
Definition: selection.h:904
void set_frame(int fr)
Set current frame for selection.
Definition: selection.cpp:601
float & X(int ind)
Extracts X for current frame.
Definition: selection.h:781
void rotate(int axis, float angle)
Rotate selection around given axis relative to center of masses.
Definition: selection.cpp:1075
bool operator!=(const Atom_proxy &other) const
Inequality operator.
Definition: selection.h:1135
float & Charge(int ind)
Extracts atom charge.
Definition: selection.h:913
void apply_transform(const Eigen::Affine3f &t)
Apply fitting transformation.
std::vector< std::string > get_resname() const
Get vector of all resnames in selection.
Definition: selection.cpp:805
bool operator==(const Atom_proxy &other) const
Equality operator.
Definition: selection.h:1130
void set_mass(const std::vector< float > m)
Set atom masses in selection to the values from supplied vector.
Definition: selection.cpp:757
bool coord_dependent() const
Returns true if selection is coordinate-dependent and is able to recompute itself on the change of fr...
Definition: selection.h:702
float distance(int i, int j, bool is_periodic=true, Vector3i_const_ref dims=Eigen::Vector3i::Ones()) const
Get distance between two atoms (periodic in given dimensions if needed).
Definition: selection.cpp:1852
void split_by_chain(std::vector< Selection > &chains)
Split selection by chain.
Definition: selection.cpp:1724
Selection operator~() const
Creates new Selection, which is a logical negation of existing one.
Definition: selection.cpp:482
std::string gromacs_ndx(std::string name)
Returns a string formatted as Gromacs ndx file containing the current selection with given name...
Definition: selection.cpp:1557
void set_chain(const std::vector< char > &data)
Set chains from supplied vector.
Definition: selection.cpp:650
void unwrap(Vector3i_const_ref dims=Eigen::Vector3i::Ones())
Unwraps selection to make it whole if possible (without jumps over periodic box boundary).
Definition: selection.cpp:1873
std::vector< float > get_mass() const
Get masses of all atoms in selection.
Definition: selection.cpp:717
void set_beta(std::vector< float > &data)
Set beta in selection to the values from supplied vector.
Definition: selection.cpp:891
void each_residue(std::vector< Selection > &sel) const
Selects each residue, which is referenced by selection.
Definition: selection.cpp:1577
std::string & Resname(int ind)
Extracts residue name.
Definition: selection.h:877
std::vector< int >::const_iterator index_begin() const
Get const iterator for begin of index.
Definition: selection.h:317
Eigen::MatrixXf get_xyz() const
Get coordinates of all atoms in this selection for the current frame.
Definition: selection.cpp:835
int & Index(int ind)
Extracts atom index in the system, which is pointed by selection.
Definition: atom_proxy.h:908
friend Eigen::Affine3f fit_transform(const Selection &sel1, const Selection &sel2)
Returns fitting transformation for two given selections of the same size.
Definition: selection.cpp:1309
virtual ~Selection()
Destructor.
Definition: selection.cpp:189
float & Time()
Returns time stamp of the frame pointed by selection The same as:
Definition: selection.h:1021
Eigen::Vector3f & XYZ(int ind, int fr)
Extracts X,Y and Z for given frame frame fr.
Definition: selection.h:850
float & Y(int ind, int fr)
Extracts Y for given frame frame fr.
Definition: selection.h:808
void set_name(const std::vector< std::string > &data)
Set atom names in selection from supplied vector.
Definition: selection.cpp:786
float & Occupancy(int ind)
Extracts occupancy field.
Definition: selection.h:931
void minmax(Vector3f_ref min, Vector3f_ref max) const
Get minimal and maximal coordinates in selection.
Definition: selection.cpp:1466
int & Type(int ind)
Extracts type.
Definition: selection.h:859
float sasa(float probe_r=0.14, float *total_volume=nullptr, std::vector< float > *area_per_atom=nullptr, std::vector< float > *volume_per_atom=nullptr) const
Get the SASA. Returns area and computes volume and per-atom values if asked.
Definition: selection.cpp:2033
int & Resindex(int ind)
Extracts resindex.
Definition: selection.h:976
std::vector< int > get_resid() const
Get vector of all resid&#39;s in selection.
Definition: selection.cpp:680
void set_xyz(MatrixXf_const_ref coord)
Set coordinates of this selection for current frame.
Definition: selection.cpp:870
void split(std::vector< Selection > &parts, F callback)
Split selection by the values returned by custom callback function.
Definition: selection.h:743
std::vector< int > get_resindex() const
Get vector of all resindexes in selection.
Definition: selection.cpp:708
Eigen::Vector3f center(bool mass_weighted=false, bool periodic=false) const
Get the center of selection.
Definition: selection.cpp:975
void apply()
Recomputes selection without re-parsing selection text.
Definition: selection.cpp:594
void clear()
Clears selection and frees memory.
Definition: selection.cpp:269
void set_occupancy(std::vector< float > &data)
Set occupancy in selection to the values from supplied vector.
Definition: selection.cpp:920
bool operator==(const Selection &other) const
Equality operator Selection are compared by their indexes.
Definition: selection.cpp:474
std::vector< int > get_unique_resid() const
Get vector of unique resid&#39;s in selection.
Definition: selection.cpp:726
std::string & Name(int ind)
Extracts atom name.
Definition: selection.h:895
float & Beta(int ind)
Extracts B-factor.
Definition: selection.h:922
float & Z(int ind)
Extracts Z for current frame.
Definition: selection.h:817
std::vector< std::string > get_name() const
Get vector of all atom names in selection.
Definition: selection.cpp:776
Class encapsulating all operations with arbitrary triclinic periodic boxes This class stores the peri...
Definition: periodic_box.h:47
void update()
Recomputes selection completely.
Definition: selection.cpp:589