Pteros  2.0
Molecular modeling library for human beings!
structure.h
1 // Copyright Maarten L. Hekkelman, Radboud University 2008-2011.
2 // Distributed under the Boost Software License, Version 1.0.
3 // (See accompanying file LICENSE_1_0.txt or copy at
4 // http://www.boost.org/LICENSE_1_0.txt)
5 
6 #pragma once
7 
8 #include "primitives-3d.h"
9 #include "align-2d.h"
10 
11 #include "pteros/core/selection.h"
12 
13 struct MAtom;
14 class MResidue;
15 class MChain;
16 class MProtein;
17 
18 // forward declaration of buffer
19 //template<typename T, uint32 N> class buffer;
20 //typedef buffer<MResidue*,100> MResidueQueue;
21 
22 const uint32 kHistogramSize = 30;
23 
24 // a limited set of known atoms. This is an obvious candidate for improvement of DSSP.
25 enum MAtomType
26 {
27  kUnknownAtom,
28  kHydrogen,
29  // ...
30  kCarbon,
31  kNitrogen,
32  kOxygen,
33  kFluorine,
34  // ...
35  kPhosphorus,
36  kSulfur,
37  kChlorine,
38  kMagnesium,
39  kPotassium,
40  kCalcium,
41  kZinc,
42  kSelenium,
43 
44  kAtomTypeCount
45 };
46 
47 MAtomType MapElement(std::string inElement);
48 
49 // for now, MAtom contains exactly what the ATOM line contains in a PDB file
50 struct MAtom
51 {
52  uint32 mSerial;
53  char mName[5];
54  char mAltLoc;
55  char mResName[5];
56  char mChainID;
57  int16 mResSeq;
58  char mICode;
59  MAtomType mType;
60  MPoint mLoc;
61  double mOccupancy;
62  double mTempFactor;
63  char mElement[3];
64  int mCharge;
65 
66  void SetChainID(char inID) { mChainID = inID;}
67  std::string GetName() const { return mName; }
68  void Translate(const MPoint& inTranslation) { mLoc += inTranslation; }
69  void Rotate(const MQuaternion& inRotation) { mLoc.Rotate(inRotation); }
70  void WritePDB(std::ostream& os) const;
71 
72  operator const MPoint&() const { return mLoc; }
73  operator MPoint&() { return mLoc; }
74 };
75 
76 enum MResidueType
77 {
78  kUnknownResidue,
79 
80  //
81  kAlanine, // A ala
82  kArginine, // R arg
83  kAsparagine, // N asn
84  kAsparticAcid, // D asp
85  kCysteine, // C cys
86  kGlutamicAcid, // E glu
87  kGlutamine, // Q gln
88  kGlycine, // G gly
89  kHistidine, // H his
90  kIsoleucine, // I ile
91  kLeucine, // L leu
92  kLysine, // K lys
93  kMethionine, // M met
94  kPhenylalanine, // F phe
95  kProline, // P pro
96  kSerine, // S ser
97  kThreonine, // T thr
98  kTryptophan, // W trp
99  kTyrosine, // Y tyr
100  kValine, // V val
101 
102  kResidueTypeCount
103 };
104 
105 struct MResidueInfo
106 {
107  MResidueType type;
108  char code;
109  char name[4];
110 };
111 
112 // a residue number to info mapping
113 extern const MResidueInfo kResidueInfo[];
114 
115 MResidueType MapResidue(std::string inName);
116 
117 struct HBond
118 {
119  MResidue* residue;
120  double energy;
121 };
122 
123 enum MBridgeType
124 {
125  btNoBridge, btParallel, btAntiParallel
126 };
127 
128 struct MBridgeParner
129 {
130  MResidue* residue;
131  uint32 ladder;
132  bool parallel;
133 };
134 
135 enum MHelixFlag
136 {
137  helixNone, helixStart, helixEnd, helixStartAndEnd, helixMiddle
138 };
139 
140 enum MSecondaryStructure
141 {
142  loop, //' '
143  alphahelix, // H
144  betabridge, // B
145  strand, // E
146  helix_3, // G
147  helix_5, // I
148  turn, // T
149  bend // S
150 };
151 
152 class MResidue
153 {
154  public:
155  MResidue(const MResidue& residue);
156  MResidue(uint32 inNumber, char inTypeCode, MResidue* inPrevious);
157  MResidue(uint32 inNumber,
158  MResidue* inPrevious, const std::vector<MAtom>& inAtoms);
159 
160  void SetChainID(char inID);
161  char GetChainID() const { return mChainID; }
162 
163  MResidueType GetType() const { return mType; }
164 
165  const MAtom& GetCAlpha() const { return mCA; }
166  const MAtom& GetC() const { return mC; }
167  const MAtom& GetN() const { return mN; }
168  const MAtom& GetO() const { return mO; }
169  const MAtom& GetH() const { return mH; }
170 
171  double Phi() const;
172  double Psi() const;
173  std::tr1::tuple<double,char>
174  Alpha() const;
175  double Kappa() const;
176  double TCO() const;
177 
178  double Accessibility() const { return mAccessibility; }
179 
180  void SetSecondaryStructure(MSecondaryStructure inSS)
181  { mSecondaryStructure = inSS; }
182  MSecondaryStructure GetSecondaryStructure() const { return mSecondaryStructure; }
183 
184  const MResidue* Next() const { return mNext; }
185  const MResidue* Prev() const { return mPrev; }
186 
187  void SetPrev(MResidue* inResidue);
188 
189  void SetBetaPartner(uint32 n, MResidue* inResidue, uint32 inLadder,
190  bool inParallel);
191  MBridgeParner GetBetaPartner(uint32 n) const;
192 
193  void SetSheet(uint32 inSheet) { mSheet = inSheet; }
194  uint32 GetSheet() const { return mSheet; }
195 
196  bool IsBend() const { return mBend; }
197  void SetBend(bool inBend) { mBend = inBend; }
198 
199  MHelixFlag GetHelixFlag(uint32 inHelixStride) const;
200  bool IsHelixStart(uint32 inHelixStride) const;
201  void SetHelixFlag(uint32 inHelixStride, MHelixFlag inHelixFlag);
202 
203  void SetSSBridgeNr(uint8 inBridgeNr);
204  uint8 GetSSBridgeNr() const;
205 
206  void AddAtom(MAtom& inAtom);
207 
208  HBond* Donor() { return mHBondDonor; }
209  HBond* Acceptor() { return mHBondAcceptor; }
210 
211  const HBond* Donor() const { return mHBondDonor; }
212  const HBond* Acceptor() const { return mHBondAcceptor; }
213 
214  bool ValidDistance(const MResidue& inNext) const;
215 
216  static bool TestBond(const MResidue* a, const MResidue* b)
217  {
218  return a->TestBond(b);
219  }
220 
221  // bridge functions
222  MBridgeType TestBridge(MResidue* inResidue) const;
223 
224  uint16 GetSeqNumber() const { return mSeqNumber; }
225  char GetInsertionCode() const { return mInsertionCode; }
226 
227  void SetNumber(uint16 inNumber) { mNumber = inNumber; }
228  uint16 GetNumber() const { return mNumber; }
229 
230  void Translate(const MPoint& inTranslation);
231  void Rotate(const MQuaternion& inRotation);
232 
233  void WritePDB(std::ostream& os);
234 
235  static double CalculateHBondEnergy(MResidue& inDonor, MResidue& inAcceptor);
236 
237  std::vector<MAtom>& GetSideChain() { return mSideChain; }
238  const std::vector<MAtom>&
239  GetSideChain() const { return mSideChain; }
240 
241  void GetPoints(std::vector<MPoint>& outPoints) const;
242 
243  void CalculateSurface(const std::vector<MResidue*>& inResidues);
244 
245  void GetCenterAndRadius(MPoint& outCenter, double& outRadius) const
246  { outCenter = mCenter; outRadius = mRadius; }
247 
248  static bool NoChainBreak(const MResidue* from, const MResidue* to);
249 
250  protected:
251 
252  double CalculateSurface(
253  const MAtom& inAtom, double inRadius,
254  const std::vector<MResidue*>& inResidues);
255 
256  bool TestBond(const MResidue* other) const;
257 
258  void ExtendBox(const MAtom& atom, double inRadius);
259  bool AtomIntersectsBox(const MAtom& atom, double inRadius) const;
260 
261  char mChainID;
262  MResidue* mPrev;
263  MResidue* mNext;
264  int32 mSeqNumber, mNumber;
265  char mInsertionCode;
266  MResidueType mType;
267  uint8 mSSBridgeNr;
268  double mAccessibility;
269  MSecondaryStructure mSecondaryStructure;
270  MAtom mC, mN, mCA, mO, mH;
271  HBond mHBondDonor[2], mHBondAcceptor[2];
272  std::vector<MAtom> mSideChain;
273  MBridgeParner mBetaPartner[2];
274  uint32 mSheet;
275  MHelixFlag mHelixFlags[3]; //
276  bool mBend;
277  MPoint mBox[2]; // The 3D box containing all atoms
278  MPoint mCenter; // and the 3d Sphere containing all atoms
279  double mRadius;
280 
281  private:
282  MResidue& operator=(const MResidue& residue);
283 };
284 
285 class MChain
286 {
287  public:
288 
289  MChain(const MChain& chain);
290  MChain(char inChainID = 0) : mChainID(inChainID) {}
291  ~MChain();
292 
293  MChain& operator=(const MChain& chain);
294 
295  char GetChainID() const { return mChainID; }
296  void SetChainID(char inID);
297 
298  MResidue* GetResidueBySeqNumber(uint16 inSeqNumber, char inInsertionCode);
299 
300  void GetSequence(std::string& outSequence) const;
301 
302  void Translate(const MPoint& inTranslation);
303  void Rotate(const MQuaternion& inRotation);
304 
305  void WritePDB(std::ostream& os);
306 
307  std::vector<MResidue*>&
308  GetResidues() { return mResidues; }
309  const std::vector<MResidue*>&
310  GetResidues() const { return mResidues; }
311 
312  bool Empty() const { return mResidues.empty(); }
313 
314  private:
315  char mChainID;
316  std::vector<MResidue*>
317  mResidues;
318 };
319 
320 class MProtein
321 {
322  public:
323  MProtein() {}
324  MProtein(const std::string& inID, MChain* inChain);
325  ~MProtein();
326 
327  const std::string& GetID() const { return mID; }
328 
329  MProtein(std::istream& is, bool inCAlphaOnly = false);
330 
331  //===========
332  MProtein(pteros::Selection& sel);
333  //===========
334 
335  const std::string& GetHeader() const { return mHeader; }
336  std::string GetCompound() const;
337  std::string GetSource() const;
338  std::string GetAuthor() const;
339  const std::vector<std::string>&
340  GetDbRef() const { return mDbRef; }
341 
342  void CalculateSecondaryStructure(bool inPreferPiHelices = true);
343 
344  void GetStatistics(uint32& outNrOfResidues, uint32& outNrOfChains,
345  uint32& outNrOfSSBridges, uint32& outNrOfIntraChainSSBridges,
346  uint32& outNrOfHBonds, uint32 outNrOfHBondsPerDistance[11]) const;
347 
348  void GetCAlphaLocations(char inChain, std::vector<MPoint>& outPoints) const;
349  MPoint GetCAlphaPosition(char inChain, int16 inPDBResSeq) const;
350 
351  void GetSequence(char inChain, entry& outEntry) const;
352  void GetSequence(char inChain, sequence& outSequence) const;
353 
354  void Center();
355  void Translate(const MPoint& inTranslation);
356  void Rotate(const MQuaternion& inRotation);
357 
358  void WritePDB(std::ostream& os);
359 
360  void GetPoints(std::vector<MPoint>& outPoints) const;
361 
362  char GetFirstChainID() const { return mChains.front()->GetChainID(); }
363 
364  void SetChain(char inChainID, const MChain& inChain);
365 
366  MChain& GetChain(char inChainID);
367  const MChain& GetChain(char inChainID) const;
368 
369  const std::vector<MChain*>&
370  GetChains() const { return mChains; }
371 
372  template<class OutputIterator>
373  void GetSequences(OutputIterator outSequences) const;
374 
375  MResidue* GetResidue(char inChainID, uint16 inSeqNumber, char inInsertionCode);
376 
377  // statistics
378  uint32 GetNrOfHBondsInParallelBridges() const { return mNrOfHBondsInParallelBridges; }
379  uint32 GetNrOfHBondsInAntiparallelBridges() const { return mNrOfHBondsInAntiparallelBridges; }
380 
381  void GetResiduesPerAlphaHelixHistogram(uint32 outHistogram[30]) const;
382  void GetParallelBridgesPerLadderHistogram(uint32 outHistogram[30]) const;
383  void GetAntiparallelBridgesPerLadderHistogram(uint32 outHistogram[30]) const;
384  void GetLaddersPerSheetHistogram(uint32 outHistogram[30]) const;
385 
386  private:
387 
388  void AddResidue(const std::vector<MAtom>& inAtoms);
389 
390  void CalculateHBondEnergies(const std::vector<MResidue*>& inResidues);
391  void CalculateAlphaHelices(const std::vector<MResidue*>& inResidues, bool inPreferPiHelices);
392  void CalculateBetaSheets(const std::vector<MResidue*>& inResidues);
393 
394  /* !! Not ised in Pteros !!
395  void CalculateAccessibilities(const std::vector<MResidue*>& inResidues);
396 
397 
398  // a thread entry point
399  void CalculateAccessibility(MResidueQueue& inQueue,
400  const std::vector<MResidue*>& inResidues);
401  */
402 
403  std::string mID, mHeader;
404 
405  std::vector<std::string>
406  mCompound, mSource, mAuthor, mDbRef;
407  std::vector<MChain*>mChains;
408  uint32 mResidueCount, mChainBreaks;
409 
410  std::vector<std::pair<MResidue*,MResidue*> >
411  mSSBonds;
412  uint32 mIgnoredWaterMolecules;
413 
414  // statistics
415  uint32 mNrOfHBondsInParallelBridges, mNrOfHBondsInAntiparallelBridges;
416  uint32 mParallelBridgesPerLadderHistogram[kHistogramSize];
417  uint32 mAntiparallelBridgesPerLadderHistogram[kHistogramSize];
418  uint32 mLaddersPerSheetHistogram[kHistogramSize];
419 };
420 
421 // inlines
422 
423 // GetSequences can be used to quickly get all sequences in a vector<string> e.g.
424 template<class OutputIterator>
425 void MProtein::GetSequences(OutputIterator outSequences) const
426 {
427  for (std::vector<MChain*>::const_iterator chain = mChains.begin(); chain != mChains.end(); ++chain)
428  {
429  std::string seq;
430  (*chain)->GetSequence(seq);
431  *outSequences++ = seq;
432  }
433 }
434 
Selection class.
Definition: atom_proxy.h:57