OpenStructure
 All Data Structures Namespaces Files Functions Variables Typedefs Enumerations Enumerator Properties Friends Macros Groups Pages
buildingblock.hh
Go to the documentation of this file.
1 //------------------------------------------------------------------------------
2 // This file is part of the OpenStructure project <www.openstructure.org>
3 //
4 // Copyright (C) 2008-2015 by the OpenStructure authors
5 //
6 // This library is free software; you can redistribute it and/or modify it under
7 // the terms of the GNU Lesser General Public License as published by the Free
8 // Software Foundation; either version 3.0 of the License, or (at your option)
9 // any later version.
10 // This library is distributed in the hope that it will be useful, but WITHOUT
11 // ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
12 // FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public License for more
13 // details.
14 //
15 // You should have received a copy of the GNU Lesser General Public License
16 // along with this library; if not, write to the Free Software Foundation, Inc.,
17 // 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
18 //------------------------------------------------------------------------------
19 
20 #ifndef OST_MM_BUILDING_BLOCK_HH
21 #define OST_MM_BUILDING_BLOCK_HH
22 
23 #include <vector>
24 #include <limits>
25 #include <set>
26 
27 #include <boost/shared_ptr.hpp>
28 
30 #include <ost/message.hh>
31 #include <ost/mol/bond_handle.hh>
33 #include <ost/mol/atom_handle.hh>
34 #include <ost/mol/xcs_editor.hh>
35 #include <ost/geom/vec3.hh>
36 
37 
38 namespace ost { namespace mol{ namespace mm{
39 
41 typedef boost::shared_ptr<BuildingBlock> BuildingBlockPtr;
42 
44 public:
46 
47  //copy constructor needs to be defined explicitely, since the forces
48  //have to be newly created
49  BuildingBlock(const BuildingBlock& block);
50 
51  bool Match(const ost::mol::ResidueHandle& handle, bool match_connectivity, String& info) const;
52 
54 
55  //getter functionality for all buildingblock members
56 
57  std::vector<String> GetAtoms() const { return atoms_; }
58 
59  std::vector<String> GetTypes() const { return types_; }
60 
61  std::vector<Real> GetCharges() const { return charges_; }
62 
63  std::vector<Real> GetMasses() const { return masses_; }
64 
65  String GetType(const String& name) const;
66 
67  Real GetCharge(const String& name) const;
68 
69  Real GetMass(const String& name) const;
70 
71  std::vector<InteractionPtr> GetBonds() const { return bonds_; }
72 
73  std::vector<InteractionPtr> GetAngles() const { return angles_; }
74 
75  std::vector<InteractionPtr> GetDihedrals() const { return dihedrals_; }
76 
77  std::vector<InteractionPtr> GetImpropers() const { return impropers_; }
78 
79  std::vector<InteractionPtr> GetCMaps() const { return cmaps_; }
80 
81  std::vector<InteractionPtr> GetExclusions() const { return exclusions_; }
82 
83  std::vector<InteractionPtr> GetConstraints() const { return constraints_;}
84 
85  //Add data to building block
86 
87  void AddAtom(const String& name, const String& type, Real charge, Real mass = std::numeric_limits<Real>::quiet_NaN());
88 
89  void AddBond(InteractionPtr p, bool replace_existing = false);
90 
91  void AddAngle(InteractionPtr p, bool replace_existing = false);
92 
93  void AddDihedral(InteractionPtr p, bool replace_existing = false);
94 
95  void AddImproper(InteractionPtr p, bool replace_existing = false);
96 
97  void AddExclusion(InteractionPtr p, bool replace_existing = false);
98 
99  void AddCMap(InteractionPtr p, bool replace_existing = false);
100 
101  void AddConstraint(InteractionPtr p, bool replace_existing = false);
102 
103  //modifiers
104 
105  //removes atom and all interactions associated to it
106  void RemoveAtom(const String& name);
107 
108  //replaces atom in all interactions
109  void ReplaceAtom(const String& name,const String& new_name,
110  const String& new_type, Real new_charge,
111  Real new_mass = std::numeric_limits<Real>::quiet_NaN());
112 
113 
114  //remove all interactions to previous or next residue
116 
118 
119  template <typename DS>
120  void Serialize(DS& ds){
121  int num_atoms = atoms_.size();
122  ds & num_atoms;
123 
124  if(ds.IsSource()){
125  atoms_ = std::vector<String>(num_atoms);
126  types_ = std::vector<String>(num_atoms);
127  charges_ = std::vector<Real>(num_atoms);
128  masses_ = std::vector<Real>(num_atoms);
129  }
130 
131  for(int i = 0; i < num_atoms; ++i){
132  ds & atoms_[i];
133  ds & types_[i];
134  ds & charges_[i];
135  ds & masses_[i];
136  }
137 
138  int num_bonds = bonds_.size();
139  int num_angles = angles_.size();
140  int num_dihedrals = dihedrals_.size();
141  int num_impropers = impropers_.size();
142  int num_exclusions = exclusions_.size();
143  int num_cmaps = cmaps_.size();
144  int num_constraints = constraints_.size();
145 
146  ds & num_bonds;
147  ds & num_angles;
148  ds & num_dihedrals;
149  ds & num_impropers;
150  ds & num_exclusions;
151  ds & num_cmaps;
152  ds & num_constraints;
153 
154  for(int i = 0; i < num_bonds; ++i){
155  int func_type;
156  if(ds.IsSource()){
157  ds & func_type;
158  bonds_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
159  }
160  else{
161  func_type = bonds_[i]->GetFuncType();
162  ds & func_type;
163  }
164  ds & *(bonds_[i]);
165  }
166 
167  for(int i = 0; i < num_angles; ++i){
168  int func_type;
169  if(ds.IsSource()){
170  ds & func_type;
171  angles_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
172  }
173  else{
174  func_type = angles_[i]->GetFuncType();
175  ds & func_type;
176  }
177  ds & *(angles_[i]);
178  }
179 
180  for(int i = 0; i < num_dihedrals; ++i){
181  int func_type;
182  if(ds.IsSource()){
183  ds & func_type;
184  dihedrals_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
185  }
186  else{
187  func_type = dihedrals_[i]->GetFuncType();
188  ds & func_type;
189  }
190  ds & *(dihedrals_[i]);
191  }
192 
193  for(int i = 0; i < num_impropers; ++i){
194  int func_type;
195  if(ds.IsSource()){
196  ds & func_type;
197  impropers_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
198  }
199  else{
200  func_type = impropers_[i]->GetFuncType();
201  ds & func_type;
202  }
203  ds & *(impropers_[i]);
204  }
205 
206  for(int i = 0; i < num_exclusions; ++i){
207  int func_type;
208  if(ds.IsSource()){
209  ds & func_type;
210  exclusions_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
211  }
212  else{
213  func_type = exclusions_[i]->GetFuncType();
214  ds & func_type;
215  }
216  ds & *(exclusions_[i]);
217  }
218 
219  for(int i = 0; i < num_cmaps; ++i){
220  int func_type;
221  if(ds.IsSource()){
222  ds & func_type;
223  cmaps_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
224  }
225  else{
226  func_type = cmaps_[i]->GetFuncType();
227  ds & func_type;
228  }
229  ds & *(cmaps_[i]);
230  }
231 
232  for(int i = 0; i < num_constraints; ++i){
233  int func_type;
234  if(ds.IsSource()){
235  ds & func_type;
236  constraints_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
237  }
238  else{
239  func_type = constraints_[i]->GetFuncType();
240  ds & func_type;
241  }
242  ds & *(constraints_[i]);
243  }
244  }
245 
246 
247 private:
248 
249  int GetAtomIndex(const String& atom_name) const;
250  void CheckInteractionToAdd(InteractionPtr p) const;
251 
252  std::vector<String> atoms_;
253  std::vector<String> types_;
254  std::vector<Real> charges_;
255  std::vector<Real> masses_;
256  std::vector<InteractionPtr> bonds_;
257  std::vector<InteractionPtr> angles_;
258  std::vector<InteractionPtr> dihedrals_;
259  std::vector<InteractionPtr> impropers_;
260  std::vector<InteractionPtr> exclusions_;
261  std::vector<InteractionPtr> cmaps_;
262  std::vector<InteractionPtr> constraints_;
263 };
264 
265 }}}
266 
267 #endif
String GetType(const String &name) const
void AddConstraint(InteractionPtr p, bool replace_existing=false)
void AddCMap(InteractionPtr p, bool replace_existing=false)
std::vector< InteractionPtr > GetDihedrals() const
void AddAtom(const String &name, const String &type, Real charge, Real mass=std::numeric_limits< Real >::quiet_NaN())
std::vector< String > GetAtoms() const
boost::shared_ptr< Interaction > InteractionPtr
Definition: interaction.hh:39
std::string String
Definition: base.hh:54
float Real
Definition: base.hh:44
std::vector< Real > GetMasses() const
Real GetCharge(const String &name) const
void AddExclusion(InteractionPtr p, bool replace_existing=false)
void ReplaceAtom(const String &name, const String &new_name, const String &new_type, Real new_charge, Real new_mass=std::numeric_limits< Real >::quiet_NaN())
void AddBond(InteractionPtr p, bool replace_existing=false)
std::vector< String > GetTypes() const
std::vector< InteractionPtr > GetAngles() const
void AddImproper(InteractionPtr p, bool replace_existing=false)
std::vector< Real > GetCharges() const
std::vector< InteractionPtr > GetBonds() const
void AddAngle(InteractionPtr p, bool replace_existing=false)
std::vector< InteractionPtr > GetImpropers() const
external coordinate system editor
Definition: xcs_editor.hh:36
Real GetMass(const String &name) const
std::vector< InteractionPtr > GetConstraints() const
bool Match(const ost::mol::ResidueHandle &handle, bool match_connectivity, String &info) const
void AddDihedral(InteractionPtr p, bool replace_existing=false)
std::vector< InteractionPtr > GetCMaps() const
std::vector< InteractionPtr > GetExclusions() const
void Connect(ost::mol::ResidueHandle &handle, ost::mol::XCSEditor &ed)
void RemoveAtom(const String &name)
boost::shared_ptr< BuildingBlock > BuildingBlockPtr