00001
00002
00003
00004
00005
00006
00007
00008
00009
00010
00011
00012
00013
00014
00015
00016
00017
00018
00019
00020 #ifndef OST_MM_BUILDING_BLOCK_HH
00021 #define OST_MM_BUILDING_BLOCK_HH
00022
00023 #include <vector>
00024 #include <limits>
00025 #include <set>
00026
00027 #include <boost/shared_ptr.hpp>
00028
00029 #include <ost/mol/mm/interaction.hh>
00030 #include <ost/message.hh>
00031 #include <ost/mol/bond_handle.hh>
00032 #include <ost/mol/residue_handle.hh>
00033 #include <ost/mol/atom_handle.hh>
00034 #include <ost/mol/xcs_editor.hh>
00035 #include <ost/geom/vec3.hh>
00036
00037
00038 namespace ost { namespace mol{ namespace mm{
00039
00040 class BuildingBlock;
00041 typedef boost::shared_ptr<BuildingBlock> BuildingBlockPtr;
00042
00043 class BuildingBlock{
00044 public:
00045 BuildingBlock() { }
00046
00047
00048
00049 BuildingBlock(const BuildingBlock& block);
00050
00051 bool Match(const ost::mol::ResidueHandle& handle, bool match_connectivity, String& info) const;
00052
00053 void Connect(ost::mol::ResidueHandle& handle, ost::mol::XCSEditor& ed);
00054
00055
00056
00057 std::vector<String> GetAtoms() const { return atoms_; }
00058
00059 std::vector<String> GetTypes() const { return types_; }
00060
00061 std::vector<Real> GetCharges() const { return charges_; }
00062
00063 std::vector<Real> GetMasses() const { return masses_; }
00064
00065 String GetType(const String& name) const;
00066
00067 Real GetCharge(const String& name) const;
00068
00069 Real GetMass(const String& name) const;
00070
00071 std::vector<InteractionPtr> GetBonds() const { return bonds_; }
00072
00073 std::vector<InteractionPtr> GetAngles() const { return angles_; }
00074
00075 std::vector<InteractionPtr> GetDihedrals() const { return dihedrals_; }
00076
00077 std::vector<InteractionPtr> GetImpropers() const { return impropers_; }
00078
00079 std::vector<InteractionPtr> GetCMaps() const { return cmaps_; }
00080
00081 std::vector<InteractionPtr> GetExclusions() const { return exclusions_; }
00082
00083 std::vector<InteractionPtr> GetConstraints() const { return constraints_;}
00084
00085
00086
00087 void AddAtom(const String& name, const String& type, Real charge, Real mass = std::numeric_limits<Real>::quiet_NaN());
00088
00089 void AddBond(InteractionPtr p, bool replace_existing = false);
00090
00091 void AddAngle(InteractionPtr p, bool replace_existing = false);
00092
00093 void AddDihedral(InteractionPtr p, bool replace_existing = false);
00094
00095 void AddImproper(InteractionPtr p, bool replace_existing = false);
00096
00097 void AddExclusion(InteractionPtr p, bool replace_existing = false);
00098
00099 void AddCMap(InteractionPtr p, bool replace_existing = false);
00100
00101 void AddConstraint(InteractionPtr p, bool replace_existing = false);
00102
00103
00104
00105
00106 void RemoveAtom(const String& name);
00107
00108
00109 void ReplaceAtom(const String& name,const String& new_name,
00110 const String& new_type, Real new_charge,
00111 Real new_mass = std::numeric_limits<Real>::quiet_NaN());
00112
00113
00114
00115 void RemoveInteractionsToPrev();
00116
00117 void RemoveInteractionsToNext();
00118
00119 template <typename DS>
00120 void Serialize(DS& ds){
00121 int num_atoms = atoms_.size();
00122 ds & num_atoms;
00123
00124 if(ds.IsSource()){
00125 atoms_ = std::vector<String>(num_atoms);
00126 types_ = std::vector<String>(num_atoms);
00127 charges_ = std::vector<Real>(num_atoms);
00128 masses_ = std::vector<Real>(num_atoms);
00129 }
00130
00131 for(int i = 0; i < num_atoms; ++i){
00132 ds & atoms_[i];
00133 ds & types_[i];
00134 ds & charges_[i];
00135 ds & masses_[i];
00136 }
00137
00138 int num_bonds = bonds_.size();
00139 int num_angles = angles_.size();
00140 int num_dihedrals = dihedrals_.size();
00141 int num_impropers = impropers_.size();
00142 int num_exclusions = exclusions_.size();
00143 int num_cmaps = cmaps_.size();
00144 int num_constraints = constraints_.size();
00145
00146 ds & num_bonds;
00147 ds & num_angles;
00148 ds & num_dihedrals;
00149 ds & num_impropers;
00150 ds & num_exclusions;
00151 ds & num_cmaps;
00152 ds & num_constraints;
00153
00154 for(int i = 0; i < num_bonds; ++i){
00155 int func_type;
00156 if(ds.IsSource()){
00157 ds & func_type;
00158 bonds_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
00159 }
00160 else{
00161 func_type = bonds_[i]->GetFuncType();
00162 ds & func_type;
00163 }
00164 ds & *(bonds_[i]);
00165 }
00166
00167 for(int i = 0; i < num_angles; ++i){
00168 int func_type;
00169 if(ds.IsSource()){
00170 ds & func_type;
00171 angles_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
00172 }
00173 else{
00174 func_type = angles_[i]->GetFuncType();
00175 ds & func_type;
00176 }
00177 ds & *(angles_[i]);
00178 }
00179
00180 for(int i = 0; i < num_dihedrals; ++i){
00181 int func_type;
00182 if(ds.IsSource()){
00183 ds & func_type;
00184 dihedrals_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
00185 }
00186 else{
00187 func_type = dihedrals_[i]->GetFuncType();
00188 ds & func_type;
00189 }
00190 ds & *(dihedrals_[i]);
00191 }
00192
00193 for(int i = 0; i < num_impropers; ++i){
00194 int func_type;
00195 if(ds.IsSource()){
00196 ds & func_type;
00197 impropers_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
00198 }
00199 else{
00200 func_type = impropers_[i]->GetFuncType();
00201 ds & func_type;
00202 }
00203 ds & *(impropers_[i]);
00204 }
00205
00206 for(int i = 0; i < num_exclusions; ++i){
00207 int func_type;
00208 if(ds.IsSource()){
00209 ds & func_type;
00210 exclusions_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
00211 }
00212 else{
00213 func_type = exclusions_[i]->GetFuncType();
00214 ds & func_type;
00215 }
00216 ds & *(exclusions_[i]);
00217 }
00218
00219 for(int i = 0; i < num_cmaps; ++i){
00220 int func_type;
00221 if(ds.IsSource()){
00222 ds & func_type;
00223 cmaps_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
00224 }
00225 else{
00226 func_type = cmaps_[i]->GetFuncType();
00227 ds & func_type;
00228 }
00229 ds & *(cmaps_[i]);
00230 }
00231
00232 for(int i = 0; i < num_constraints; ++i){
00233 int func_type;
00234 if(ds.IsSource()){
00235 ds & func_type;
00236 constraints_.push_back(InteractionPtr(new Interaction(FuncType(func_type))));
00237 }
00238 else{
00239 func_type = constraints_[i]->GetFuncType();
00240 ds & func_type;
00241 }
00242 ds & *(constraints_[i]);
00243 }
00244 }
00245
00246
00247 private:
00248
00249 int GetAtomIndex(const String& atom_name) const;
00250 void CheckInteractionToAdd(InteractionPtr p) const;
00251
00252 std::vector<String> atoms_;
00253 std::vector<String> types_;
00254 std::vector<Real> charges_;
00255 std::vector<Real> masses_;
00256 std::vector<InteractionPtr> bonds_;
00257 std::vector<InteractionPtr> angles_;
00258 std::vector<InteractionPtr> dihedrals_;
00259 std::vector<InteractionPtr> impropers_;
00260 std::vector<InteractionPtr> exclusions_;
00261 std::vector<InteractionPtr> cmaps_;
00262 std::vector<InteractionPtr> constraints_;
00263 };
00264
00265 }}}
00266
00267 #endif