1
/**********************************************************************
2
mol.h - Handle molecules. Declarations of OBMol, OBAtom, OBBond, OBResidue.
3
(the main header for Open Babel)
5
Copyright (C) 1998-2001 by OpenEye Scientific Software, Inc.
6
Some portions Copyright (C) 2001-2005 by Geoffrey R. Hutchison
7
Some portions Copyright (C) 2003 by Michael Banck
9
This file is part of the Open Babel project.
10
For more information, see <http://openbabel.sourceforge.net/>
12
This program is free software; you can redistribute it and/or modify
13
it under the terms of the GNU General Public License as published by
14
the Free Software Foundation version 2 of the License.
16
This program is distributed in the hope that it will be useful,
17
but WITHOUT ANY WARRANTY; without even the implied warranty of
18
MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the
19
GNU General Public License for more details.
20
***********************************************************************/
25
#include "babelconfig.h"
28
# define EXTERN extern
53
#include "math/vector3.h"
60
#include "reaction.h" //so it gets notices in DLL builds
68
class OBInternalCoord;
71
// class introduction in residue.cpp
79
//! \warning Currently does not copy all associated OBGenericData
80
//! This requires a (minor) API change, and will thus only be fixed in 2.1
81
//! or later releases.
82
OBResidue(const OBResidue &);
84
virtual ~OBResidue(void);
86
OBResidue &operator=(const OBResidue &);
88
void AddAtom(OBAtom *atom);
89
void InsertAtom(OBAtom *atom);
90
void RemoveAtom(OBAtom *atom);
93
void SetName(const std::string &resname);
94
void SetNum(unsigned int resnum);
95
void SetChain(char chain);
96
void SetChainNum(unsigned int chainnum);
97
void SetIdx(unsigned int idx);
99
void SetAtomID(OBAtom *atom, const std::string &id);
100
void SetHetAtom(OBAtom *atom, bool hetatm);
101
//! Set the atomic serial number for a given atom (see OBSerialNums)
102
void SetSerialNum(OBAtom *atom, unsigned int sernum);
104
std::string GetName(void) const;
105
unsigned int GetNum(void) const;
106
unsigned int GetNumAtoms() const;
107
char GetChain(void) const;
108
unsigned int GetChainNum(void) const;
109
unsigned int GetIdx(void) const;
110
unsigned int GetResKey(void) const;
112
std::vector<OBAtom*> GetAtoms(void) const;
113
std::vector<OBBond*> GetBonds(bool = true) const;
115
std::string GetAtomID(OBAtom *atom) const;
116
//! \return the serial number of the supplied atom (uses OBSerialNums)
117
unsigned GetSerialNum(OBAtom *atom) const;
119
bool GetAminoAcidProperty(int) const;
120
bool GetAtomProperty(OBAtom *, int) const;
121
bool GetResidueProperty(int) const;
123
bool IsHetAtom(OBAtom *atom) const;
124
bool IsResidueType(int) const;
126
//! \deprecated Use FOR_ATOMS_OF_RESIDUE and OBResidueAtomIter instead
127
OBAtom *BeginAtom(std::vector<OBAtom*>::iterator &i);
128
//! \deprecated Use FOR_ATOMS_OF_RESIDUE and OBResidueAtomIter instead
129
OBAtom *NextAtom(std::vector<OBAtom*>::iterator &i);
131
//! \name Methods for handling generic data
133
bool HasData(std::string &);
134
bool HasData(const char *);
135
bool HasData(unsigned int type);
136
void DeleteData(unsigned int type);
137
void DeleteData(OBGenericData*);
138
void DeleteData(std::vector<OBGenericData*>&);
139
void SetData(OBGenericData *d)
140
{ _vdata.push_back(d); }
141
//! \return the number of OBGenericData items attached to this residue
142
unsigned int DataSize()
143
{ return(_vdata.size()); }
144
OBGenericData *GetData(unsigned int type);
145
OBGenericData *GetData(std::string&);
146
OBGenericData *GetData(const char *);
147
std::vector<OBGenericData*> &GetData()
149
std::vector<OBGenericData*>::iterator BeginData()
150
{ return(_vdata.begin()); }
151
std::vector<OBGenericData*>::iterator EndData()
152
{ return(_vdata.end()); }
155
protected: // members
157
unsigned int _idx; //!< Residue index (i.e., internal index in an OBMol)
158
char _chain; //!< Chain ID
159
unsigned int _aakey; //!< Amino Acid key ID -- see SetResidueKeys()
160
unsigned int _reskey;//!< Residue key ID -- see SetResidueKeys()
161
unsigned int _resnum;//!< Residue number (i.e., in file)
162
std::string _resname;//!< Residue text name
164
std::vector<bool> _hetatm;//!< Is a given atom a HETAM
165
std::vector<std::string> _atomid;//!< Residue atom text IDs
166
std::vector<OBAtom*> _atoms; //!< List of OBAtom in this residue
167
std::vector<unsigned int> _sernum;//!< List of serial numbers
168
std::vector<OBGenericData*> _vdata; //!< Custom data
172
//ATOM Property Macros (flags)
173
//! Atom is in a 4-membered ring
174
#define OB_4RING_ATOM (1<<1)
175
//! Atom is in a 3-membered ring
176
#define OB_3RING_ATOM (1<<2)
178
#define OB_AROMATIC_ATOM (1<<3)
179
//! Atom is in a ring
180
#define OB_RING_ATOM (1<<4)
181
//! Atom has clockwise SMILES chiral stereochemistry (i.e., "@@")
182
#define OB_CSTEREO_ATOM (1<<5)
183
//! Atom has anticlockwise SMILES chiral stereochemistry (i.e., "@")
184
#define OB_ACSTEREO_ATOM (1<<6)
185
//! Atom is an electron donor
186
#define OB_DONOR_ATOM (1<<7)
187
//! Atom is an electron acceptor
188
#define OB_ACCEPTOR_ATOM (1<<8)
190
#define OB_CHIRAL_ATOM (1<<9)
191
//! Atom has + chiral volume
192
#define OB_POS_CHIRAL_ATOM (1<<10)
193
//! Atom has - chiral volume
194
#define OB_NEG_CHIRAL_ATOM (1<<11)
195
// 12-16 currently unused
198
// class introduction in atom.cpp
199
class OBAPI OBAtom : public OBNodeBase
202
char _ele; //!< atomic number (type char to minimize space -- allows for 0..255 elements)
203
char _impval; //!< implicit valence
204
char _type[6]; //!< atomic type
205
short _fcharge; //!< formal charge
206
unsigned short _isotope; //!< isotope (0 = most abundant)
207
short _spinmultiplicity;//!< atomic spin, e.g., 2 for radical 1 or 3 for carbene
209
//unsigned short int _idx; //!< index in parent (inherited)
210
unsigned short _cidx; //!< index into coordinate array
211
unsigned short _hyb; //!< hybridization
212
unsigned short _flags; //!< bitwise flags (e.g. aromaticity)
213
double _pcharge; //!< partial charge
214
double **_c; //!< coordinate array in double*
215
vector3 _v; //!< coordinate vector
216
OBResidue *_residue; //!< parent residue (if applicable)
217
//OBMol *_parent; //!< parent molecule (inherited)
218
//vector<OBBond*> _bond; //!< connections (inherited)
219
std::vector<OBGenericData*> _vdata; //!< custom data
221
int GetFlag() const { return(_flags); }
222
void SetFlag(int flag) { _flags |= flag; }
223
bool HasFlag(int flag) { return((_flags & flag) ? true : false); }
232
OBAtom &operator = (OBAtom &);
236
//! \name Methods to set atomic information
238
//! Set atom index (i.e., in an OBMol)
239
void SetIdx(int idx) { _idx = idx; _cidx = (idx-1)*3; }
240
//! Set atom hybridization (i.e., 1 = sp, 2 = sp2, 3 = sp3 ...)
241
void SetHyb(int hyb) { _hyb = hyb; }
242
//! Set atomic number
243
void SetAtomicNum(int atomicnum) { _ele = (char)atomicnum; }
244
//! Set isotope number (actual atomic weight is tabulated automatically, 0 = most abundant)
245
void SetIsotope(unsigned int iso);
246
void SetImplicitValence(int val) { _impval = (char)val; }
247
void IncrementImplicitValence() { _impval++; }
248
void DecrementImplicitValence() { _impval--; }
249
void SetFormalCharge(int fcharge) { _fcharge = fcharge; }
250
void SetSpinMultiplicity(short spin){ _spinmultiplicity = spin; }
251
void SetType(char *type);
252
void SetType(std::string &type);
253
void SetPartialCharge(double pcharge){ _pcharge = pcharge; }
254
void SetVector(vector3 &v);
255
void SetVector(const double x,const double y,const double z);
256
//! Set the position of this atom from a pointer-driven array of coordinates
257
void SetCoordPtr(double **c) { _c = c; _cidx = (GetIdx()-1)*3; }
258
//! Set the position of this atom based on the internal pointer array (i.e. from SetCoordPtr() )
260
void SetResidue(OBResidue *res) { _residue=res; }
261
// void SetParent(OBMol *ptr) { _parent=ptr; } // inherited
262
void SetAromatic() { SetFlag(OB_AROMATIC_ATOM); }
263
void UnsetAromatic() { _flags &= (~(OB_AROMATIC_ATOM)); }
264
//! Mark atom as having SMILES clockwise stereochemistry (i.e., "@@")
265
void SetClockwiseStereo() { SetFlag(OB_CSTEREO_ATOM|OB_CHIRAL_ATOM); }
266
//! Mark atom as having SMILES anticlockwise stereochemistry (i.e., "@")
267
void SetAntiClockwiseStereo() { SetFlag(OB_ACSTEREO_ATOM|OB_CHIRAL_ATOM); }
268
//! Mark an atom as having + chiral volume
269
void SetPositiveStereo() { SetFlag(OB_POS_CHIRAL_ATOM|OB_CHIRAL_ATOM); }
270
//! Mark an atom as having - chiral volume
271
void SetNegativeStereo() { SetFlag(OB_NEG_CHIRAL_ATOM|OB_CHIRAL_ATOM); }
272
//! Clear all stereochemistry information
275
_flags &= ~(OB_ACSTEREO_ATOM);
276
_flags &= ~(OB_CSTEREO_ATOM);
277
_flags &= ~(OB_POS_CHIRAL_ATOM);
278
_flags &= ~(OB_NEG_CHIRAL_ATOM);
279
_flags &= ~(OB_CHIRAL_ATOM);
281
//! Mark an atom as belonging to at least one ring
282
void SetInRing() { SetFlag(OB_RING_ATOM); }
283
//! Mark an atom as being chiral with unknown stereochemistry
284
void SetChiral() { SetFlag(OB_CHIRAL_ATOM); }
285
//! Clear the internal coordinate pointer
286
void ClearCoordPtr() { _c = NULL; _cidx=0; }
289
//! \name Methods to retrieve atomic information
291
//int GetStereo() const { return((int)_stereo);}
292
int GetFormalCharge() const { return(_fcharge); }
293
unsigned int GetAtomicNum() const { return((unsigned int)_ele); }
294
unsigned short int GetIsotope() const { return(_isotope); }
295
int GetSpinMultiplicity() const { return(_spinmultiplicity); }
296
//! The atomic mass of this atom given by standard IUPAC average molar mass
297
double GetAtomicMass() const;
298
//! The atomic mass of given by the isotope (default of 0 s most abundant isotope)
299
double GetExactMass() const;
300
unsigned int GetIdx() const { return((int)_idx); }
301
unsigned int GetCoordinateIdx() const { return((int)_cidx); }
302
//! \deprecated Use GetCoordinateIdx() instead
303
unsigned int GetCIdx() const { return((int)_cidx); }
304
//! The current number of explicit connections
305
unsigned int GetValence() const
307
return((_vbond.empty()) ? 0 : _vbond.size());
309
//! The hybridization of this atom (i.e. 1 for sp, 2 for sp2, 3 for sp3)
310
unsigned int GetHyb() const;
311
//! The implicit valence of this atom type (i.e. maximum number of connections expected)
312
unsigned int GetImplicitValence() const;
313
//! The number of non-hydrogens connected to this atom
314
unsigned int GetHvyValence() const;
315
//! The number of heteroatoms connected to an atom
316
unsigned int GetHeteroValence() const;
320
double GetX() { return(x()); }
322
double GetY() { return(y()); }
324
double GetZ() { return(z()); }
328
return((*_c)[_cidx]);
335
return((*_c)[_cidx+1]);
342
return((*_c)[_cidx+2]);
346
//! \return the coordinates as a double*
347
double *GetCoordinate()
350
return(&(*_c)[_cidx]);
354
//! \return the coordinates as a vector3 object
355
vector3 &GetVector();
356
//! \return the partial charge of this atom, calculating a Gasteiger charge if needed
357
double GetPartialCharge();
358
OBResidue *GetResidue();
359
//OBMol *GetParent() {return((OBMol*)_parent);}
360
//! Create a vector for a new bond from this atom, with length given by the supplied parameter
361
bool GetNewBondVector(vector3 &v,double length);
362
OBBond *GetBond(OBAtom *);
363
OBAtom *GetNextAtom();
366
//! \name Iterator methods
368
//! \deprecated Use FOR_BONDS_OF_ATOM and OBAtomBondIter instead
369
std::vector<OBEdgeBase*>::iterator BeginBonds()
370
{ return(_vbond.begin()); }
371
//! \deprecated Use FOR_BONDS_OF_ATOM and OBAtomBondIter instead
372
std::vector<OBEdgeBase*>::iterator EndBonds()
373
{ return(_vbond.end()); }
374
//! \deprecated Use FOR_BONDS_OF_ATOM and OBAtomBondIter instead
375
OBBond *BeginBond(std::vector<OBEdgeBase*>::iterator &i);
376
//! \deprecated Use FOR_BONDS_OF_ATOM and OBAtomBondIter instead
377
OBBond *NextBond(std::vector<OBEdgeBase*>::iterator &i);
378
//! \deprecated Use FOR_NBORS_OF_ATOM and OBAtomAtomIter instead
379
OBAtom *BeginNbrAtom(std::vector<OBEdgeBase*>::iterator &);
380
//! \deprecated Use FOR_NBORS_OF_ATOM and OBAtomAtomIter instead
381
OBAtom *NextNbrAtom(std::vector<OBEdgeBase*>::iterator &);
384
//! \return the distance to the atom defined by OBMol::GetAtom()
385
double GetDistance(int index);
386
//! \return the distance to the supplied OBAtom
387
double GetDistance(OBAtom*);
388
//! \return the angle defined by this atom -> b (vertex) -> c
389
double GetAngle(int b, int c);
390
//! \return the angle defined by this atom -> b (vertex) -> c
391
double GetAngle(OBAtom *b, OBAtom *c);
393
//! \name Addition of residue/bond info. for an atom
398
_residue = new OBResidue;
405
void AddBond(OBBond *bond)
407
_vbond.push_back((OBEdgeBase*)bond);
409
void InsertBond(std::vector<OBEdgeBase*>::iterator &i, OBBond *bond)
411
_vbond.insert(i, (OBEdgeBase*)bond);
413
bool DeleteBond(OBBond*);
414
void ClearBond() {_vbond.clear();}
417
//! \name Requests for atomic property information
419
//! The number of oxygen atoms connected that only have one heavy valence
420
unsigned int CountFreeOxygens() const;
421
//! The number of hydrogens needed to fill the implicit valence of this atom
422
unsigned int ImplicitHydrogenCount() const;
423
//! The number of hydrogens explicitly bound to this atom currently
424
unsigned int ExplicitHydrogenCount() const;
425
//! The number of rings that contain this atom
426
unsigned int MemberOfRingCount() const;
427
//! The size of the smallest ring that contains this atom (0 if not in a ring)
428
unsigned int MemberOfRingSize() const;
429
//! The smallest angle of bonds to this atom
430
double SmallestBondAngle();
431
//! The average angle of bonds to this atom
432
double AverageBondAngle();
433
//! The sum of the bond orders of the bonds to the atom (i.e. double bond = 2...)
434
unsigned int BOSum() const;
435
//! The sum of the bond orders of bonds to the atom, considering only KDouble, KTriple bonds
436
unsigned int KBOSum() const;
439
//! \name Builder utilities
441
//! If this is a hydrogen atom, transform into a methyl group
443
//! Change the hybridization of this atom and modify the geometry accordingly
444
bool SetHybAndGeom(int);
447
//! \name Property information
449
//! Is there any residue information?
450
bool HasResidue() { return(_residue != NULL); }
451
bool IsHydrogen() { return(GetAtomicNum() == 1); }
452
bool IsCarbon() { return(GetAtomicNum() == 6); }
453
bool IsNitrogen() { return(GetAtomicNum() == 7); }
454
bool IsOxygen() { return(GetAtomicNum() == 8); }
455
bool IsSulfur() { return(GetAtomicNum() == 16);}
456
bool IsPhosphorus() { return(GetAtomicNum() == 15);}
457
bool IsAromatic() const;
458
bool IsInRing() const;
459
bool IsInRingSize(int) const;
460
//! Is this atom an element in the 15th or 16th main groups (i.e., N, O, P, S ...) ?
462
//! Is this atom any element except carbon or hydrogen?
464
//! Is this atom connected to the supplied OBAtom?
465
bool IsConnected(OBAtom*);
466
//! Is this atom related to the supplied OBAtom in a 1,3 bonding pattern?
467
bool IsOneThree(OBAtom*);
468
//! Is this atom related to the supplied OBAtom in a 1,4 bonding pattern?
469
bool IsOneFour(OBAtom*);
470
//! Is this atom an oxygen in a carboxyl (-CO2 or CO2H) group?
471
bool IsCarboxylOxygen();
472
//! Is this atom an oxygen in a phosphate (R-PO3) group?
473
bool IsPhosphateOxygen();
474
//! Is this atom an oxygen in a sulfate (-SO3) group?
475
bool IsSulfateOxygen();
476
//! Is this atom an oxygen in a nitro (-NO2) group?
477
bool IsNitroOxygen();
478
bool IsAmideNitrogen();
479
bool IsPolarHydrogen();
480
bool IsNonPolarHydrogen();
481
bool IsAromaticNOxide();
482
//! Is this atom chiral?
485
//! Does this atom have SMILES-specified clockwise "@@" stereochemistry?
486
bool IsClockwise() { return(HasFlag(OB_CSTEREO_ATOM)); }
487
//! Does this atom have SMILES-specified anticlockwise "@" stereochemistry?
488
bool IsAntiClockwise() { return(HasFlag(OB_ACSTEREO_ATOM)); }
489
//! Does this atom have a positive chiral volume?
490
bool IsPositiveStereo() { return(HasFlag(OB_POS_CHIRAL_ATOM)); }
491
//! Does this atom have a negative chiral volume?
492
bool IsNegativeStereo() { return(HasFlag(OB_NEG_CHIRAL_ATOM)); }
493
//! Does this atom have SMILES-specified stereochemistry?
494
bool HasChiralitySpecified()
495
{ return(HasFlag(OB_CSTEREO_ATOM|OB_ACSTEREO_ATOM)); }
496
//! Does this atom have a specified chiral volume?
497
bool HasChiralVolume()
498
{ return(HasFlag(OB_POS_CHIRAL_ATOM|OB_NEG_CHIRAL_ATOM)); }
499
//! Is this atom a hydrogen-bond acceptor (receptor)?
500
bool IsHbondAcceptor();
501
//! Is this atom a hydrogen-bond donor?
503
//! Is this a hydrogen atom attached to a hydrogen-bond donor?
504
bool IsHbondDonorH();
505
bool HasAlphaBetaUnsat(bool includePandS=true);
506
bool HasBondOfOrder(unsigned int);
507
int CountBondsOfOrder(unsigned int);
508
bool HasNonSingleBond();
509
bool HasSingleBond() { return(HasBondOfOrder(1)); }
510
bool HasDoubleBond() { return(HasBondOfOrder(2)); }
511
bool HasAromaticBond() { return(HasBondOfOrder(5)); }
512
//! Determines if this atom matches the first atom in a given SMARTS pattern
513
bool MatchesSMARTS(const char *);
516
//! \name Methods for handling generic data
518
bool HasData(std::string &);
519
bool HasData(const char *);
520
bool HasData(unsigned int type);
521
void DeleteData(unsigned int type);
522
void DeleteData(OBGenericData*);
523
void DeleteData(std::vector<OBGenericData*>&);
524
void SetData(OBGenericData *d)
525
{ _vdata.push_back(d); }
526
//! \return the number of OBGenericData items attached to this atom
527
unsigned int DataSize()
528
{ return(_vdata.size()); }
529
OBGenericData *GetData(unsigned int type);
530
OBGenericData *GetData(std::string&);
531
OBGenericData *GetData(const char *);
532
std::vector<OBGenericData*> &GetData() { return(_vdata); }
533
std::vector<OBGenericData*>::iterator BeginData()
534
{ return(_vdata.begin()); }
535
std::vector<OBGenericData*>::iterator EndData()
536
{ return(_vdata.end()); }
543
//BOND Property Macros (flags)
544
//! An aromatic bond (regardless of bond order)
545
#define OB_AROMATIC_BOND (1<<1)
546
//! A solid black wedge in 2D representations -- i.e., "up" from the 2D plane
547
#define OB_WEDGE_BOND (1<<2)
548
//! A dashed "hash" bond in 2D representations -- i.e., "down" from the 2D plane
549
#define OB_HASH_BOND (1<<3)
551
#define OB_RING_BOND (1<<4)
552
//! The "upper" bond in a double bond cis/trans isomer (i.e., "/" in SMILES)
553
#define OB_TORUP_BOND (1<<5)
554
//! The "down" bond in a double bond cis/trans isomer (i.e., "\" in SMILES)
555
#define OB_TORDOWN_BOND (1<<6)
556
//! A Kekule single bond
557
#define OB_KSINGLE_BOND (1<<7)
558
//! A Kekule double bond
559
#define OB_KDOUBLE_BOND (1<<8)
560
//! A Kekule triple bond
561
#define OB_KTRIPLE_BOND (1<<9)
562
#define OB_CLOSURE_BOND (1<<10)
563
// 11-16 currently unused
565
// class introduction in bond.cpp
566
class OBAPI OBBond : public OBEdgeBase
569
char _order; //!< Bond order (1, 2, 3, 5=aromatic)
570
unsigned short int _flags; //!< Any flags for this bond
571
//OBAtom *_bgn; //!< Not needed, inherited from OBEdgeBase
572
//OBAtom *_end; //!< Not needed, inherited from OBEdgeBase
573
//OBMol *_parent;//!< Not needed, inherited from OBEdgeBase
574
//unsigned short int _idx; //!< Not needed, inherited from OBEdgeBase
575
std::vector<OBGenericData*> _vdata; //!< Generic data for custom information
577
bool HasFlag(int flag) { return((_flags & flag) != 0); }
578
void SetFlag(int flag) { _flags |= flag; }
586
//! \name Bond modification methods
592
void SetBO(int order);
593
void SetBegin(OBAtom *begin)
597
void SetEnd(OBAtom *end)
601
// void SetParent(OBMol *ptr) {_parent=ptr;} // (inherited)
602
void SetLength(OBAtom*,double);
603
void Set(int,OBAtom*,OBAtom*,int,int);
607
void SetAromatic() { SetFlag(OB_AROMATIC_BOND); }
608
void SetHash() { SetFlag(OB_HASH_BOND); }
609
void SetWedge() { SetFlag(OB_WEDGE_BOND); }
610
void SetUp() { SetFlag(OB_TORUP_BOND); }
611
void SetDown() { SetFlag(OB_TORDOWN_BOND); }
612
void SetInRing() { SetFlag(OB_RING_BOND); }
613
void SetClosure() { SetFlag(OB_CLOSURE_BOND); }
615
void UnsetAromatic() { _flags &= (~(OB_AROMATIC_BOND)); }
618
_flags &= (~(OB_KSINGLE_BOND|OB_KDOUBLE_BOND|OB_KTRIPLE_BOND));
622
//! \name bond data request methods
624
unsigned int GetBO() const { return((int)_order); }
625
unsigned int GetBondOrder() const { return((int)_order); }
626
unsigned int GetFlags() const { return(_flags); }
627
unsigned int GetBeginAtomIdx() const { return(_bgn->GetIdx()); }
628
unsigned int GetEndAtomIdx() const { return(_end->GetIdx()); }
629
OBAtom *GetBeginAtom() { return((OBAtom*)_bgn); }
630
OBAtom *GetEndAtom() { return((OBAtom*)_end); }
631
OBAtom *GetNbrAtom(OBAtom *ptr)
633
return((ptr != _bgn)? (OBAtom*)_bgn : (OBAtom*)_end);
635
// OBMol *GetParent() {return(_parent);} // (inherited)
636
double GetEquibLength();
638
int GetNbrAtomIdx(OBAtom *ptr)
640
return((ptr!=_bgn)?_bgn->GetIdx():_end->GetIdx());
644
//! \name property request methods
646
bool IsAromatic() const;
647
bool IsInRing() const;
648
//! Is the bond a rotatable bond?
649
//! Currently, this function classifies any bond with at least one heavy
650
//! atom, no sp-hybrid atoms (e.g., a triple bond somewhere) not in a ring
651
//! as a potential rotor. No other bond typing is attempted.
654
bool IsPrimaryAmide();
655
bool IsSecondaryAmide();
665
//! \return whether this is the "upper" bond in a double bond cis/trans
666
//! isomer (i.e., "/" in SMILES)
667
bool IsUp() { return(HasFlag(OB_TORUP_BOND)); }
668
//! \return whether this is the "lower" bond in a double bond cis/trans
669
//! isomer (i.e., "\" in SMILES)
670
bool IsDown() { return(HasFlag(OB_TORDOWN_BOND)); }
671
bool IsWedge() { return(HasFlag(OB_WEDGE_BOND)); }
672
bool IsHash() { return(HasFlag(OB_HASH_BOND)); }
673
//! \return whether the geometry around this bond looks unsaturated
674
bool IsDoubleBondGeometry();
677
//! \name Methods for handling generic data
679
bool HasData(std::string &);
680
bool HasData(const char *);
681
bool HasData(unsigned int type);
682
void DeleteData(unsigned int type);
683
void DeleteData(OBGenericData*);
684
void DeleteData(std::vector<OBGenericData*>&);
685
void SetData(OBGenericData *d)
689
//! \return the number of OBGenericData items attached to this bond
690
unsigned int DataSize()
692
return(_vdata.size());
694
OBGenericData *GetData(unsigned int type);
695
OBGenericData *GetData(std::string&);
696
OBGenericData *GetData(const char *);
697
std::vector<OBGenericData*> &GetData()
701
std::vector<OBGenericData*>::iterator BeginData()
703
return(_vdata.begin());
705
std::vector<OBGenericData*>::iterator EndData()
707
return(_vdata.end());
716
//MOL Property Macros (flags) -- 32+ bits
717
#define OB_SSSR_MOL (1<<1)
718
#define OB_RINGFLAGS_MOL (1<<2)
719
#define OB_AROMATIC_MOL (1<<3)
720
#define OB_ATOMTYPES_MOL (1<<4)
721
#define OB_CHIRALITY_MOL (1<<5)
722
#define OB_PCHARGE_MOL (1<<6)
723
#define OB_HYBRID_MOL (1<<8)
724
#define OB_IMPVAL_MOL (1<<9)
725
#define OB_KEKULE_MOL (1<<10)
726
#define OB_CLOSURE_MOL (1<<11)
727
#define OB_H_ADDED_MOL (1<<12)
728
#define OB_PH_CORRECTED_MOL (1<<13)
729
#define OB_AROM_CORRECTED_MOL (1<<14)
730
#define OB_CHAINS_MOL (1<<15)
731
#define OB_TCHARGE_MOL (1<<16)
732
#define OB_TSPIN_MOL (1<<17)
733
// flags 18-32 unspecified
734
#define OB_CURRENT_CONFORMER -1
736
// class introduction in mol.cpp
737
class OBAPI OBMol : public OBGraphBase
740
int _flags; //!< bitfield of flags
741
bool _autoPartialCharge; //!< Assign partial charges automatically
742
bool _autoFormalCharge; //!< Assign formal charges automatically
743
std::string _title; //!< Molecule title
744
//vector<OBAtom*> _atom; //!< not needed (inherited)
745
//vector<OBBond*> _bond; //!< not needed (inherited)
746
unsigned short int _dimension; //!< Dimensionality of coordinates
747
double _energy; //!< Molecular heat of formation (if applicable)
748
int _totalCharge; //!< Total charge on the molecule
749
unsigned int _totalSpin; //!< Total spin on the molecule (if not specified, assumes lowest possible spin)
750
double *_c; //!< coordinate array
751
std::vector<double*> _vconf; //!< vector of conformers
752
unsigned short int _natoms; //!< Number of atoms
753
unsigned short int _nbonds; //!< Number of bonds
754
std::vector<OBResidue*> _residue; //!< Residue information (if applicable)
755
std::vector<OBInternalCoord*> _internals; //!< Internal Coordinates (if applicable)
756
std::vector<OBGenericData*> _vdata; //!< Custom data -- see OBGenericData class for more
757
unsigned short int _mod; //!< Number of nested calls to BeginModify()
759
bool HasFlag(int flag) { return((_flags & flag) ? true : false); }
760
void SetFlag(int flag) { _flags |= flag; }
762
//! \name Internal Kekulization routines -- see kekulize.cpp and NewPerceiveKekuleBonds()
764
void start_kekulize(std::vector <OBAtom*> &cycle, std::vector<int> &electron);
765
int expand_kekulize(OBAtom *atom1, OBAtom *atom2, std::vector<int> ¤tState, std::vector<int> &initState, std::vector<int> &bcurrentState, std::vector<int> &binitState, std::vector<bool> &mark);
766
int getorden(OBAtom *atom);
767
void expandcycle(OBAtom *atom, OBBitVec &avisit);
772
//! \name Initialization and data (re)size methods
777
//! \warning Currently does not copy all associated OBGenericData
778
//! This requires a (minor) API change, and will thus only be fixed in 2.1
779
//! or later releases.
780
OBMol(const OBMol &);
784
OBMol &operator=(const OBMol &mol);
785
OBMol &operator+=(const OBMol &mol);
786
void ReserveAtoms(int natoms)
789
_vatom.reserve(natoms);
791
virtual OBAtom *CreateAtom(void);
792
virtual OBBond *CreateBond(void);
793
virtual void DestroyAtom(OBNodeBase*);
794
virtual void DestroyBond(OBEdgeBase*);
795
bool AddAtom(OBAtom&);
796
bool AddBond(int,int,int,int flags=0,int insertpos=-1);
797
bool AddBond(OBBond&);
798
bool AddResidue(OBResidue&);
799
bool InsertAtom(OBAtom &);
800
bool DeleteAtom(OBAtom*);
801
bool DeleteBond(OBBond*);
802
bool DeleteResidue(OBResidue*);
804
OBResidue *NewResidue();
807
//! \name Molecule modification methods
809
//! Call when making many modifications -- clears conformer/rotomer data.
810
virtual void BeginModify(void);
811
//! Call when done with modificaions -- re-perceive data as needed.
812
virtual void EndModify(bool nukePerceivedData=true);
827
//! \name Generic data handling methods (via OBGenericData)
829
//! \returns whether the generic attribute/value pair exists
830
bool HasData(std::string &);
831
//! \returns whether the generic attribute/value pair exists
832
bool HasData(const char *);
833
//! \returns whether the generic attribute/value pair exists
834
bool HasData(unsigned int type);
835
void DeleteData(unsigned int type);
836
void DeleteData(OBGenericData*);
837
void DeleteData(std::vector<OBGenericData*>&);
838
void SetData(OBGenericData *d)
842
//! \return the number of OBGenericData items attached to this molecule.
843
unsigned int DataSize(){ return(_vdata.size()); }
844
OBGenericData *GetData(unsigned int type);
845
OBGenericData *GetData(std::string&);
846
OBGenericData *GetData(const char *);
847
std::vector<OBGenericData*> &GetData() { return(_vdata); }
848
std::vector<OBGenericData*>::iterator BeginData()
850
return(_vdata.begin());
852
std::vector<OBGenericData*>::iterator EndData()
854
return(_vdata.end());
858
//! \name Data retrieval methods
860
int GetFlags() { return(_flags); }
861
//! \return the title of this molecule (often the filename)
862
const char *GetTitle() const { return(_title.c_str()); }
863
//! \return the number of atoms (i.e. OBAtom children)
864
unsigned int NumAtoms() const { return(_natoms); }
865
//! \return the number of bonds (i.e. OBBond children)
866
unsigned int NumBonds() const { return(_nbonds); }
867
//! \return the number of non-hydrogen atoms
868
unsigned int NumHvyAtoms();
869
//! \return the number of residues (i.e. OBResidue substituents)
870
unsigned int NumResidues() const { return(_residue.size()); }
871
//! \return the number of rotatble bonds. See OBBond::IsRotor() for details
872
unsigned int NumRotors();
874
OBAtom *GetAtom(int);
875
OBAtom *GetFirstAtom();
876
OBBond *GetBond(int);
877
OBBond *GetBond(int, int);
878
OBBond *GetBond(OBAtom*,OBAtom*);
879
OBResidue *GetResidue(int);
880
std::vector<OBInternalCoord*> GetInternalCoord();
881
//! \return the dihedral angle between the four atoms supplied a1-a2-a3-a4)
882
double GetTorsion(int,int,int,int);
883
//! \return the dihedral angle between the four atoms supplied a1-a2-a3-a4)
884
double GetTorsion(OBAtom*,OBAtom*,OBAtom*,OBAtom*);
885
//! \return the stochoimetric formula (e.g., C4H6O)
886
std::string GetFormula();
887
//! \return the heat of formation for this molecule (in kcal/mol)
888
double GetEnergy() const { return(_energy); }
889
//! \return the standard molar mass given by IUPAC atomic masses (amu)
891
//! \return the mass given by isotopes (or most abundant isotope, if not specified)
892
double GetExactMass();
893
//! \return the total charge on this molecule (i.e., 0 = neutral, +1, -1...)
894
int GetTotalCharge();
895
//! \return the total spin on this molecule (i.e., 1 = singlet, 2 = doublet...)
896
unsigned int GetTotalSpinMultiplicity();
897
//! \return the dimensionality of coordinates (i.e., 0 = unknown or no coord, 2=2D, 3=3D)
898
unsigned short int GetDimension() const { return _dimension; }
899
double *GetCoordinates() { return(_c); }
900
//! \return the Smallest Set of Smallest Rings has been run (see OBRing class
901
std::vector<OBRing*> &GetSSSR();
902
//! Get the current flag for whether formal charges are set with pH correction
903
bool AutomaticFormalCharge() { return(_autoFormalCharge); }
904
//! Get the current flag for whether partial charges are auto-determined
905
bool AutomaticPartialCharge() { return(_autoPartialCharge); }
909
//! \name Data modification methods
911
void SetTitle(const char *title);
912
void SetTitle(std::string &title);
913
//! Set the stochiometric formula for this molecule
914
void SetFormula(std::string molFormula);
915
//! Set the heat of formation for this molecule (in kcal/mol)
916
void SetEnergy(double energy) { _energy = energy; }
917
//! Set the dimension of this molecule (i.e., 0, 1 , 2, 3)
918
void SetDimension(unsigned short int d) { _dimension = d; }
919
void SetTotalCharge(int charge);
920
void SetTotalSpinMultiplicity(unsigned int spin);
921
void SetInternalCoord(std::vector<OBInternalCoord*> int_coord)
922
{ _internals = int_coord; }
923
//! Set the flag for determining automatic formal charges with pH (default=true)
924
void SetAutomaticFormalCharge(bool val)
925
{ _autoFormalCharge=val; }
926
//! Set the flag for determining partial charges automatically (default=true)
927
void SetAutomaticPartialCharge(bool val)
928
{ _autoPartialCharge=val; }
930
//! Mark that aromaticity has been perceived for this molecule (see OBAromaticTyper)
931
void SetAromaticPerceived() { SetFlag(OB_AROMATIC_MOL); }
932
//! Mark that Smallest Set of Smallest Rings has been run (see OBRing class)
933
void SetSSSRPerceived() { SetFlag(OB_SSSR_MOL); }
934
//! Mark that rings have been perceived (see OBRing class for details)
935
void SetRingAtomsAndBondsPerceived(){SetFlag(OB_RINGFLAGS_MOL);}
936
//! Mark that atom types have been perceived (see OBAtomTyper for details)
937
void SetAtomTypesPerceived() { SetFlag(OB_ATOMTYPES_MOL); }
938
//! Mark that chains and residues have been perceived (see OBChainsParser)
939
void SetChainsPerceived() { SetFlag(OB_CHAINS_MOL); }
940
//! Mark that chirality has been perceived
941
void SetChiralityPerceived() { SetFlag(OB_CHIRALITY_MOL); }
942
//! Mark that partial charges have been assigned
943
void SetPartialChargesPerceived(){ SetFlag(OB_PCHARGE_MOL); }
944
void SetHybridizationPerceived() { SetFlag(OB_HYBRID_MOL); }
945
void SetImplicitValencePerceived(){ SetFlag(OB_IMPVAL_MOL); }
946
void SetKekulePerceived() { SetFlag(OB_KEKULE_MOL); }
947
void SetClosureBondsPerceived(){ SetFlag(OB_CLOSURE_MOL); }
948
void SetHydrogensAdded() { SetFlag(OB_H_ADDED_MOL); }
949
void SetCorrectedForPH() { SetFlag(OB_PH_CORRECTED_MOL);}
950
void SetAromaticCorrected() { SetFlag(OB_AROM_CORRECTED_MOL);}
951
void SetSpinMultiplicityAssigned(){ SetFlag(OB_TSPIN_MOL); }
952
void SetFlags(int flags) { _flags = flags; }
954
void UnsetAromaticPerceived() { _flags &= (~(OB_AROMATIC_MOL)); }
955
void UnsetPartialChargesPerceived(){ _flags &= (~(OB_PCHARGE_MOL));}
956
void UnsetImplicitValencePerceived(){_flags &= (~(OB_IMPVAL_MOL)); }
957
void UnsetFlag(int flag) { _flags &= (~(flag)); }
959
//! \name Molecule modification methods
961
// Description in transform.cpp
962
virtual OBBase* DoTransformations(const std::map<std::string,std::string>* pOptions);
963
static const char* ClassDescription();
964
//! Clear all information from a molecule
966
//! Renumber the atoms of this molecule according to the order in the supplied vector
967
void RenumberAtoms(std::vector<OBNodeBase*>&);
968
//! Translate one conformer and rotate by a rotation matrix (which is returned) to the inertial frame-of-reference
969
void ToInertialFrame(int conf, double *rmat);
970
//! Translate all conformers to the inertial frame-of-reference
971
void ToInertialFrame();
972
//! Translates all conformers in the molecule by the supplied vector
973
void Translate(const vector3 &v);
974
//! Translates one conformer in the molecule by the supplied vector
975
void Translate(const vector3 &v, int conf);
976
void Rotate(const double u[3][3]);
977
void Rotate(const double m[9]);
978
void Rotate(const double m[9],int nconf);
979
//! Translate to the center of all coordinates (for this conformer)
981
//! Transform to standard Kekule bond structure (presumably from an aromatic form)
983
bool PerceiveKekuleBonds();
985
void NewPerceiveKekuleBonds();
987
bool DeleteHydrogen(OBAtom*);
988
bool DeleteHydrogens();
989
bool DeleteHydrogens(OBAtom*);
990
bool DeleteNonPolarHydrogens();
991
bool AddHydrogens(bool polaronly=false,bool correctForPH=true);
992
bool AddHydrogens(OBAtom*);
993
bool AddPolarHydrogens();
995
//! Deletes all atoms except for the largest contiguous fragment
997
//! Converts the charged form of coordinate bonds, e.g.[N+]([O-])=O to N(=O)=O
998
bool ConvertDativeBonds();
1000
bool CorrectForPH();
1001
bool AssignSpinMultiplicity();
1002
vector3 Center(int nconf);
1003
//! Set the torsion defined by these atoms, rotating bonded neighbors
1004
void SetTorsion(OBAtom*,OBAtom*,OBAtom*,OBAtom*,double);
1007
//! \name Molecule utilities and perception methods
1009
//! Find Smallest Set of Smallest Rings (see OBRing class for more details)
1011
void FindRingAtomsAndBonds();
1012
void FindChiralCenters();
1013
void FindChildren(std::vector<int> &,int,int);
1014
void FindChildren(std::vector<OBAtom*>&,OBAtom*,OBAtom*);
1015
void FindLargestFragment(OBBitVec &);
1016
//! Sort a list of contig fragments by size from largest to smallest
1017
//! Each vector<int> contains the atom numbers of a contig fragment
1018
void ContigFragList(std::vector<std::vector<int> >&);
1019
//! Aligns atom a on p1 and atom b along p1->p2 vector
1020
void Align(OBAtom*,OBAtom*,vector3&,vector3&);
1021
//! Adds single bonds based on atom proximity
1022
void ConnectTheDots();
1023
//! Attempts to perceive multiple bonds based on geometries
1024
void PerceiveBondOrders();
1025
void FindTorsions();
1026
// documented in mol.cpp: graph-theoretical distance for each atom
1027
bool GetGTDVector(std::vector<int> &);
1028
// documented in mol.cpp: graph-invariant index for each atom
1029
void GetGIVector(std::vector<unsigned int> &);
1030
// documented in mol.cpp: calculate symmetry-unique identifiers
1031
void GetGIDVector(std::vector<unsigned int> &);
1034
//! \name Methods to check for existence of properties
1036
//! Are there non-zero coordinates in two dimensions (i.e. X and Y)?
1038
//! Are there non-zero coordinates in all three dimensions (i.e. X, Y, Z)?
1040
//! Are there any non-zero coordinates?
1041
bool HasNonZeroCoords();
1042
bool HasAromaticPerceived() { return(HasFlag(OB_AROMATIC_MOL)); }
1043
bool HasSSSRPerceived() { return(HasFlag(OB_SSSR_MOL)); }
1044
bool HasRingAtomsAndBondsPerceived(){return(HasFlag(OB_RINGFLAGS_MOL));}
1045
bool HasAtomTypesPerceived() { return(HasFlag(OB_ATOMTYPES_MOL));}
1046
bool HasChiralityPerceived() { return(HasFlag(OB_CHIRALITY_MOL));}
1047
bool HasPartialChargesPerceived() { return(HasFlag(OB_PCHARGE_MOL));}
1048
bool HasHybridizationPerceived() { return(HasFlag(OB_HYBRID_MOL)); }
1049
bool HasImplicitValencePerceived() { return(HasFlag(OB_IMPVAL_MOL));}
1050
bool HasKekulePerceived() { return(HasFlag(OB_KEKULE_MOL)); }
1051
bool HasClosureBondsPerceived() { return(HasFlag(OB_CLOSURE_MOL)); }
1052
bool HasChainsPerceived() { return(HasFlag(OB_CHAINS_MOL)); }
1053
bool HasHydrogensAdded() { return(HasFlag(OB_H_ADDED_MOL)); }
1054
bool HasAromaticCorrected() { return(HasFlag(OB_AROM_CORRECTED_MOL));}
1055
bool IsCorrectedForPH() { return(HasFlag(OB_PH_CORRECTED_MOL)); }
1056
bool HasSpinMultiplicityAssigned() { return(HasFlag(OB_TSPIN_MOL)); }
1057
//! Is this molecule chiral?
1059
//! Are there any atoms in this molecule?
1060
bool Empty() { return(_natoms == 0); }
1063
//! \name Multiple conformer member functions
1065
int NumConformers() { return((_vconf.empty())?0:_vconf.size()); }
1066
void SetConformers(std::vector<double*> &v);
1067
void AddConformer(double *f) { _vconf.push_back(f); }
1068
void SetConformer(int i) { _c = _vconf[i]; }
1069
void CopyConformer(double*,int);
1070
void DeleteConformer(int);
1071
double *GetConformer(int i) { return(_vconf[i]); }
1072
double *BeginConformer(std::vector<double*>::iterator&i)
1073
{ i = _vconf.begin();
1074
return((i == _vconf.end()) ? NULL:*i); }
1075
double *NextConformer(std::vector<double*>::iterator&i)
1077
return((i == _vconf.end()) ? NULL:*i); }
1078
std::vector<double*> &GetConformers() { return(_vconf); }
1081
//! \name Iterator methods
1083
//! \deprecated Use FOR_ATOMS_OF_MOL and OBMolAtomIter instead
1084
OBAtom *BeginAtom(std::vector<OBNodeBase*>::iterator &i);
1085
//! \deprecated Use FOR_ATOMS_OF_MOL and OBMolAtomIter instead
1086
OBAtom *NextAtom(std::vector<OBNodeBase*>::iterator &i);
1087
//! \deprecated Use FOR_BONDS_OF_MOL and OBMolBondIter instead
1088
OBBond *BeginBond(std::vector<OBEdgeBase*>::iterator &i);
1089
//! \deprecated Use FOR_BONDS_OF_MOL and OBMolBondIter instead
1090
OBBond *NextBond(std::vector<OBEdgeBase*>::iterator &i);
1091
//! \deprecated Use FOR_RESIDUES_OF_MOL and OBResidueIter instead
1092
OBResidue *BeginResidue(std::vector<OBResidue*>::iterator &i)
1094
i = _residue.begin();
1095
return((i == _residue.end()) ? NULL:*i);
1097
//! \deprecated Use FOR_RESIDUES_OF_MOL and OBResidueIter instead
1098
OBResidue *NextResidue(std::vector<OBResidue*>::iterator &i)
1101
return((i == _residue.end()) ? NULL:*i);
1103
OBInternalCoord *BeginInternalCoord(std::vector<OBInternalCoord*>::iterator &i)
1105
i = _internals.begin();
1106
return((i == _internals.end()) ? NULL:*i);
1108
OBInternalCoord *NextInternalCoord(std::vector<OBInternalCoord*>::iterator &i)
1111
return((i == _internals.end()) ? NULL:*i);
1115
// Removed with OBConversion framework -- see OBConversion class instead
1116
//! \name Convenience functions for I/O
1118
// friend std::ostream& operator<< ( std::ostream&, OBMol& ) ;
1119
// friend std::istream& operator>> ( std::istream&, OBMol& ) ;
1123
//! \brief Used to transform from z-matrix to cartesian coordinates.
1124
class OBAPI OBInternalCoord
1129
double _dst,_ang,_tor;
1131
OBInternalCoord(OBAtom *a=(OBAtom*)NULL,
1132
OBAtom *b=(OBAtom*)NULL,
1133
OBAtom *c=(OBAtom*)NULL)
1138
_dst = _ang = _tor = 0.0;
1142
//function prototypes
1144
OBAPI bool tokenize(std::vector<std::string>&, const char *buf, const char *delimstr=" \t\n");
1145
OBAPI bool tokenize(std::vector<std::string>&, std::string&, const char *delimstr=" \t\n", int limit=-1);
1146
//! remove leading and trailing whitespace from a string
1147
OBAPI void Trim(std::string& txt);
1148
//! \deprecated -- use OBMessageHandler class instead
1149
OBAPI void ThrowError(char *str);
1150
//! \deprecated -- use OBMessageHandler class instead
1151
OBAPI void ThrowError(std::string &str);
1152
OBAPI void CartesianToInternal(std::vector<OBInternalCoord*>&,OBMol&);
1153
OBAPI void InternalToCartesian(std::vector<OBInternalCoord*>&,OBMol&);
1154
OBAPI std::string NewExtension(std::string&,char*);
1155
// Now handled by OBConversion class
1156
// OBAPI bool SetInputType(OBMol&,std::string&);
1157
// OBAPI bool SetOutputType(OBMol&,std::string&);
1159
//global definitions
1160
//! Global OBElementTable for element properties
1161
EXTERN OBElementTable etab;
1162
//! Global OBTypeTable for translating between different atom types
1163
//! (e.g., Sybyl <-> MM2)
1164
EXTERN OBTypeTable ttab;
1165
//! Global OBIsotopeTable for isotope properties
1166
EXTERN OBIsotopeTable isotab;
1167
//! Global OBAromaticTyper for detecting aromatic atoms and bonds
1168
EXTERN OBAromaticTyper aromtyper;
1169
//! Global OBAtomTyper for marking internal valence, hybridization,
1170
//! and atom types (for internal and external use)
1171
EXTERN OBAtomTyper atomtyper;
1172
//! Global OBChainsParser for detecting macromolecular chains and residues
1173
EXTERN OBChainsParser chainsparser;
1174
//! Global OBMessageHandler error handler
1175
EXTERN OBMessageHandler obErrorLog;
1176
//! Global OBResidueData biomolecule residue database
1177
EXTERN OBResidueData resdat;
1182
#define BUFF_SIZE 32768
1186
#define EQ(a,b) (!strcmp((a), (b)))
1190
#define EQn(a,b,n) (!strncmp((a), (b), (n)))
1194
#define SQUARE(x) ((x)*(x))
1198
#define IsUnsatType(x) (EQ(x,"Car") || EQ(x,"C2") || EQ(x,"Sox") || EQ(x,"Sac") || EQ(x,"Pac") || EQ(x,"So2"))
1204
OBAPI void get_rmat(double*,double*,double*,int);
1205
OBAPI void ob_make_rmat(double mat[3][3],double rmat[9]);
1206
OBAPI void qtrfit (double *r,double *f,int size,double u[3][3]);
1207
OBAPI double superimpose(double*,double*,int);
1210
OBAPI void get_rmat(double*,double*,double*,int);
1211
OBAPI void ob_make_rmat(double mat[3][3],double rmat[9]);
1212
OBAPI void qtrfit (double *r,double *f,int size,double u[3][3]);
1213
OBAPI double superimpose(double*,double*,int);
1216
} // end namespace OpenBabel
1221
//! \brief Handle molecules. Declarations of OBMol, OBAtom, OBBond, OBResidue.
1222
//! (the main header for Open Babel)