1
/* minimol.cpp: atomic model types */
2
//C Copyright (C) 2000-2006 Kevin Cowtan and University of York
4
//L This library is free software and is distributed under the terms
5
//L and conditions of version 2.1 of the GNU Lesser General Public
6
//L Licence (LGPL) with the following additional clause:
8
//L `You may also combine or link a "work that uses the Library" to
9
//L produce a work containing portions of the Library, and distribute
10
//L that work under terms of your choice, provided that you give
11
//L prominent notice with each copy of the work that the specified
12
//L version of the Library is used in it, and that you include or
13
//L provide public access to the complete corresponding
14
//L machine-readable source code for the Library including whatever
15
//L changes were used in the work. (i.e. If you make changes to the
16
//L Library you must distribute those, but you do not need to
17
//L distribute source or object code to those portions of the work
18
//L not covered by this licence.)'
20
//L Note that this clause grants an additional right and does not impose
21
//L any additional restriction, and so does not affect compatibility
22
//L with the GNU General Public Licence (GPL). If you wish to negotiate
23
//L other terms, please contact the maintainer.
25
//L You can redistribute it and/or modify the library under the terms of
26
//L the GNU Lesser General Public License as published by the Free Software
27
//L Foundation; either version 2.1 of the License, or (at your option) any
30
//L This library is distributed in the hope that it will be useful, but
31
//L WITHOUT ANY WARRANTY; without even the implied warranty of
32
//L MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the GNU
33
//L Lesser General Public License for more details.
35
//L You should have received a copy of the CCP4 licence and/or GNU
36
//L Lesser General Public License along with this library; if not, write
37
//L to the CCP4 Secretary, Daresbury Laboratory, Warrington WA4 4AD, UK.
38
//L The GNU Lesser General Public can also be obtained by writing to the
39
//L Free Software Foundation, Inc., 59 Temple Place, Suite 330, Boston,
48
Message_ctor message_ctor_mmodel( " [MModel: constructed]" );
53
MAtom::MAtom( const clipper::Atom& atom )
55
set_element( atom.element() );
56
set_coord_orth( atom.coord_orth() );
57
set_occupancy( atom.occupancy() );
58
set_u_iso( atom.u_iso() );
59
set_u_aniso_orth( atom.u_aniso_orth() );
62
void MAtom::set_id( const String& s ) { id_ = id_tidy( s ); }
64
void MAtom::set_name( const String s, const String altconf )
66
if ( altconf != "" ) set_id( s + ":" + altconf );
70
String MAtom::id_tidy( const String& id )
72
int pos = id.find( ":" );
73
if ( pos == String::npos ) pos = id.length();
74
String name( id.substr( 0, pos ) );
75
String altc( id.substr( pos ) );
76
if ( name.length() < 4 ) {
78
if ( islower( name[1] ) )
79
name[1] = toupper( name[1] );
83
return name.substr(0,4) + altc;
86
/*! copy from other atom. mode can be MM::COPY_M, COPY_P, COPY_MP,
87
COPY_C, COPY_MC, COPY_PC, COPY_MPC, where M means copy members, P
88
means copy PropertyMananger properties, and C means copy
89
children. Children are copied with the same option. The values
90
'MEMBERS', 'PROPERTIES', 'CHILDREN' can also be used. */
91
MAtom& MAtom::copy( const MAtom& other, const MM::COPY& mode )
93
if ( mode & MM::COPY_M ) atom() = other.atom();
94
if ( mode & MM::COPY_M ) id_ = other.id_;
95
if ( mode & MM::COPY_P ) PropertyManager::copy( other );
99
bool MAtom::id_match( const String& id1, const String& id2, const MM::MODE& mode )
100
{ if ( mode == MM::UNIQUE ) return ( id1 == id2 );
101
else return ( id1.substr(0,4) == id2.substr(0,4) ); }
106
void MMonomer::set_id( const String& s )
111
void MMonomer::set_type( const String& s )
114
void MMonomer::set_seqnum( const int s, const String inscode )
116
if ( inscode != "" ) set_id( String( s, 4 ) + ":" + inscode );
117
else set_id( String( s, 4 ) );
120
Atom_list MMonomer::atom_list() const
123
for ( int a = 0; a < children.size(); a++ )
124
list.push_back( Atom( children[a] ) );
128
void MMonomer::transform( const RTop_orth rt )
129
{ for ( int i = 0; i < children.size(); i++ ) children[i].transform( rt ); }
131
/*! Creates a copy of this monomer containing only the atoms described
132
by the selection string. '*' copies all atoms.
134
The atom selection must contain an atom ID or a comma separated list
135
of atom IDs, or '*' to select all atom. Atom IDs are described in
138
The selection string must contain an atom ID or a comma separated
139
list of atom IDs. Atom IDs are described in s_mm_atom_id.
141
\param sel The selection string.
142
\param mode MM::UNIQUE forces an exact match, including alternate conformation code. MM::ANY matches every atom with the right name, ignoring alternate conformation codes.
143
\return The selection as a new monomer. */
144
MMonomer MMonomer::select( const String& sel, const MM::MODE mode ) const
146
std::vector<String> path = sel.split( "/" );
147
while ( path.size() < 1 ) path.push_back( "*" );
149
result.copy( *this, MM::COPY_MP );
150
if ( path[0].trim() == "*" ) {
151
for ( int i = 0; i < children.size(); i++ ) result.insert( children[i] );
154
std::vector<String> list = path[0].split( "," );
155
for ( int j = 0; j < list.size(); j++ ) {
156
String sid = CHILDTYPE::id_tidy( list[j] );
157
for ( int i = 0; i < children.size(); i++ )
158
if ( CHILDTYPE::id_match( sid, children[i].id(), mode ) )
159
result.insert( children[i] );
165
/*! Lookup atom by ID. If mode=UNIQUE, the alternate conformation code
166
must match, otherwise the first atom with the same name is returned.
167
\param n The atom ID. \param mode The search mode.
169
const MAtom& MMonomer::find( const String& n, const MM::MODE mode ) const
171
int i = lookup( n, mode );
172
if ( i < 0 ) Message::message(Message_fatal("MMonomer: no such atom"));
176
/*! See MMonomer::find() */
177
MAtom& MMonomer::find( const String& n, const MM::MODE mode )
179
int i = lookup( n, mode );
180
if ( i < 0 ) Message::message(Message_fatal("MMonomer: no such atom"));
184
int MMonomer::lookup( const String& str, const MM::MODE& mode ) const
186
String sid = CHILDTYPE::id_tidy( str );
187
for ( int i = 0; i < children.size(); i++ )
188
if ( CHILDTYPE::id_match( sid, children[i].id(), mode ) ) return i;
192
void MMonomer::insert( const MAtom& add, int pos )
194
if ( pos < 0 ) children.push_back( add );
195
else children.insert( children.begin() + pos, add );
198
/*! \return atoms which are in both monomers. */
199
MMonomer operator& ( const MMonomer& m1, const MMonomer& m2 )
202
result.copy( m1, MM::COPY_MP );
203
for ( int i1 = 0; i1 < m1.size(); i1++ )
204
for ( int i2 = 0; i2 < m2.size(); i2++ )
205
if ( m1[i1].id() == m2[i2].id() ) {
206
result.insert( m1[i1] );
212
/*! \return atoms which are in either monomer. */
213
MMonomer operator| ( const MMonomer& m1, const MMonomer& m2 )
216
result.copy( m1, MM::COPY_MP );
218
for ( i = 0; i < m1.size(); i++ ) {
219
for ( j = 0; j < result.size(); j++ )
220
if ( m1[i].id() == result[j].id() ) break;
221
if ( j == result.size() )
222
result.insert( m1[i] );
224
for ( i = 0; i < m2.size(); i++ ) {
225
for ( j = 0; j < result.size(); j++ )
226
if ( m2[i].id() == result[j].id() ) break;
227
if ( j == result.size() )
228
result.insert( m2[i] );
233
/*! copy from other atom. mode can be MM::COPY_M, COPY_P, COPY_MP,
234
COPY_C, COPY_MC, COPY_PC, COPY_MPC, where M means copy members, P
235
means copy PropertyMananger properties, and C means copy
236
children. Children are copied with the same option. The values
237
'MEMBERS', 'PROPERTIES', 'CHILDREN' can also be used. */
238
MMonomer& MMonomer::copy( const MMonomer& other, const MM::COPY& mode )
240
if ( mode & MM::COPY_M ) id_ = other.id_;
241
if ( mode & MM::COPY_M ) type_ = other.type_;
242
if ( mode & MM::COPY_P ) PropertyManager::copy( other );
243
if ( mode & MM::COPY_C ) {
244
children.resize( other.size() );
245
for ( int i = 0; i < size(); i++ ) children[i].copy( other[i], mode );
250
String MMonomer::id_tidy( const String& id )
252
int pos = id.find( ":" );
253
if ( pos == String::npos )
254
return String( id.i(), 4 );
256
return String( id.i(), 4 ) + id.substr( pos );
259
bool MMonomer::id_match( const String& id1, const String& id2, const MM::MODE& mode )
260
{ if ( mode == MM::UNIQUE ) return ( id1 == id2 );
261
else return ( id1.substr(0,4) == id2.substr(0,4) ); }
266
void MPolymer::set_id( const String& s ) { id_ = id_tidy( s ); }
268
Atom_list MPolymer::atom_list() const
271
for ( int m = 0; m < children.size(); m++ )
272
for ( int a = 0; a < children[m].size(); a++ )
273
list.push_back( Atom( children[m][a] ) );
277
void MPolymer::transform( const RTop_orth rt )
278
{ for ( int i = 0; i < children.size(); i++ ) children[i].transform( rt ); }
280
/*! Creates a copy of this polymer containing only the monomers and
281
atoms described by the selection string.
283
The selection string must be of the form 'X/Y' where X is a monomer
284
selection and Y is an atom selection, described under
285
MAtom::select(). The monomer selection must contain a monomer ID or
286
a comma separated list of monomer IDs, or '*' to select all
287
monomers. Monomer IDs are described in s_mm_monomer_id.
289
\param sel The selection string.
290
\param mode MM::UNIQUE forces an exact match, including insertion code.
291
MM::ANY matches any monomer with the right sequence number,
292
ignoring insertion code.
293
\return The selection as a new polymer. */
294
MPolymer MPolymer::select( const String& sel, const MM::MODE mode ) const
296
std::vector<String> path = sel.split( "/" );
297
while ( path.size() < 2 ) path.push_back( "*" );
299
result.copy( *this, MM::COPY_MP );
300
if ( path[0].trim() == "*" ) {
301
for ( int i = 0; i < children.size(); i++ )
302
result.insert( children[i].select( path[1], mode ) );
305
std::vector<String> list = path[0].split( "," );
306
for ( int j = 0; j < list.size(); j++ ) {
307
String sid = CHILDTYPE::id_tidy( list[j] );
308
for ( int i = 0; i < children.size(); i++ )
309
if ( CHILDTYPE::id_match( sid, children[i].id(), mode ) )
310
result.insert( children[i].select( path[1], mode ) );
316
/*! Lookup monomer by ID. If mode=UNIQUE, the insertion code must match,
317
otherwise the first monomer with the same sequence number is returned.
318
\param n The monomer ID. \param mode The search mode.
319
\return The monomer. */
320
const MMonomer& MPolymer::find( const String& n, const MM::MODE mode ) const
322
int i = lookup( n, mode );
323
if ( i < 0 ) Message::message(Message_fatal("MPolymer: no such monomer"));
327
/*! See MPolymer::find() */
328
MMonomer& MPolymer::find( const String& n, const MM::MODE mode )
330
int i = lookup( n, mode );
331
if ( i < 0 ) Message::message(Message_fatal("MPolymer: no such monomer"));
335
int MPolymer::lookup( const String& str, const MM::MODE& mode ) const
337
String sid = CHILDTYPE::id_tidy( str );
338
for ( int i = 0; i < children.size(); i++ )
339
if ( CHILDTYPE::id_match( sid, children[i].id(), mode ) ) return i;
343
void MPolymer::insert( const MMonomer& add, int pos )
345
if ( pos < 0 ) children.push_back( add );
346
else children.insert( children.begin() + pos, add );
349
/*! \return monomers and atoms which are in both polymers. */
350
MPolymer operator& ( const MPolymer& m1, const MPolymer& m2 )
353
result.copy( m1, MM::COPY_MP );
354
for ( int i1 = 0; i1 < m1.size(); i1++ )
355
for ( int i2 = 0; i2 < m2.size(); i2++ )
356
if ( m1[i1].id() == m2[i2].id() ) {
357
result.insert( m1[i1] & m2[i2] );
363
/*! \return monomers and atoms which are in either polymer. */
364
MPolymer operator| ( const MPolymer& m1, const MPolymer& m2 )
367
result.copy( m1, MM::COPY_MP );
369
for ( i = 0; i < m1.size(); i++ ) {
370
for ( j = 0; j < result.size(); j++ )
371
if ( m1[i].id() == result[j].id() ) break;
372
if ( j == result.size() )
373
result.insert( m1[i] );
375
result[j] = result[j] | m1[i];
377
for ( i = 0; i < m2.size(); i++ ) {
378
for ( j = 0; j < result.size(); j++ )
379
if ( m2[i].id() == result[j].id() ) break;
380
if ( j == result.size() )
381
result.insert( m2[i] );
383
result[j] = result[j] | m2[i];
388
/*! copy from other atom. mode can be MM::COPY_M, COPY_P, COPY_MP,
389
COPY_C, COPY_MC, COPY_PC, COPY_MPC, where M means copy members, P
390
means copy PropertyMananger properties, and C means copy
391
children. Children are copied with the same option. The values
392
'MEMBERS', 'PROPERTIES', 'CHILDREN' can also be used. */
393
MPolymer& MPolymer::copy( const MPolymer& other, const MM::COPY& mode )
395
if ( mode & MM::COPY_M ) id_ = other.id_;
396
if ( mode & MM::COPY_P ) PropertyManager::copy( other );
397
if ( mode & MM::COPY_C ) {
398
children.resize( other.size() );
399
for ( int i = 0; i < size(); i++ ) children[i].copy( other[i], mode );
404
String MPolymer::id_tidy( const String& id ) { return id; }
405
bool MPolymer::id_match( const String& id1, const String& id2, const MM::MODE& mode ) { return ( id1 == id2 ); }
410
Atom_list MModel::atom_list() const
413
for ( int p = 0; p < children.size(); p++ )
414
for ( int m = 0; m < children[p].size(); m++ )
415
for ( int a = 0; a < children[p][m].size(); a++ )
416
list.push_back( Atom( children[p][m][a] ) );
420
void MModel::transform( const RTop_orth rt )
421
{ for ( int i = 0; i < children.size(); i++ ) children[i].transform( rt ); }
423
/*! Creates a copy of this model containing only the polymers,
424
monomers and atoms described by the selection string.
426
The selection string must be of the form 'X/Y/Z' where X is a
427
polymer selection, Y is a monomer selection described under
428
MMonomer::select(), and Z is an atom selection described under
429
MAtom::select(). The polymer selection must contain a polymer ID or
430
a comma separated list of polymer IDs, or '*' to select all
431
polymers. Polymer IDs are described in s_mm_monomer_id.
433
See s_mm_selections for examples.
435
\param sel The selection string.
436
\param mode No effect.
437
\return The selection as a new model. */
438
MModel MModel::select( const String& sel, const MM::MODE mode ) const
440
std::vector<String> path = sel.split( "/" );
441
while ( path.size() < 3 ) path.push_back( "*" );
443
result.copy( *this, MM::COPY_MP );
444
if ( path[0].trim() == "*" ) {
445
for ( int i = 0; i < children.size(); i++ )
446
result.insert( children[i].select( path[1]+"/"+path[2], mode ) );
449
std::vector<String> list = path[0].split( "," );
450
for ( int j = 0; j < list.size(); j++ ) {
451
String sid = CHILDTYPE::id_tidy( list[j] );
452
for ( int i = 0; i < children.size(); i++ )
453
if ( CHILDTYPE::id_match( sid, children[i].id(), mode ) )
454
result.insert( children[i].select( path[1]+"/"+path[2], mode ) );
460
/*! Lookup polymer by ID. Currently, mode is ignored.
461
\param n The monomer ID. \param mode The search mode.
462
\return The polymer. */
463
const MPolymer& MModel::find( const String& n, const MM::MODE mode ) const
465
int i = lookup( n, mode );
466
if ( i < 0 ) Message::message(Message_fatal("MModel: no such polymer"));
470
/*! See MModel::find() */
471
MPolymer& MModel::find( const String& n, const MM::MODE mode )
473
int i = lookup( n, mode );
474
if ( i < 0 ) Message::message(Message_fatal("MModel: no such polymer"));
478
int MModel::lookup( const String& str, const MM::MODE& mode ) const
480
String sid = CHILDTYPE::id_tidy( str );
481
for ( int i = 0; i < children.size(); i++ )
482
if ( CHILDTYPE::id_match( sid, children[i].id(), mode ) ) return i;
486
void MModel::insert( const MPolymer& add, int pos )
488
if ( pos < 0 ) children.push_back( add );
489
else children.insert( children.begin() + pos, add );
492
/*! \return polymers, monomers and atoms which are in both models. */
493
MModel operator& ( const MModel& m1, const MModel& m2 )
496
result.copy( m1, MM::COPY_MP );
497
for ( int i1 = 0; i1 < m1.size(); i1++ )
498
for ( int i2 = 0; i2 < m2.size(); i2++ )
499
if ( m1[i1].id() == m2[i2].id() ) {
500
result.insert( m1[i1] & m2[i2] );
506
/*! \return polymers, monomers and atoms which are in either model. */
507
MModel operator| ( const MModel& m1, const MModel& m2 )
510
result.copy( m1, MM::COPY_MP );
512
for ( i = 0; i < m1.size(); i++ ) {
513
for ( j = 0; j < result.size(); j++ )
514
if ( m1[i].id() == result[j].id() ) break;
515
if ( j == result.size() )
516
result.insert( m1[i] );
518
result[j] = result[j] | m1[i];
520
for ( i = 0; i < m2.size(); i++ ) {
521
for ( j = 0; j < result.size(); j++ )
522
if ( m2[i].id() == result[j].id() ) break;
523
if ( j == result.size() )
524
result.insert( m2[i] );
526
result[j] = result[j] | m2[i];
531
/*! copy from other atom. mode can be MM::COPY_M, COPY_P, COPY_MP,
532
COPY_C, COPY_MC, COPY_PC, COPY_MPC, where M means copy members, P
533
means copy PropertyMananger properties, and C means copy
534
children. Children are copied with the same option. The values
535
'MEMBERS', 'PROPERTIES', 'CHILDREN' can also be used. */
536
MModel& MModel::copy( const MModel& other, const MM::COPY& mode )
538
if ( mode & MM::COPY_P ) PropertyManager::copy( other );
539
if ( mode & MM::COPY_C ) {
540
children.resize( other.size() );
541
for ( int i = 0; i < size(); i++ ) children[i].copy( other[i], mode );
550
{ Message::message( message_ctor_mmodel ); }
552
/*! The object is constructed with no atoms.
553
\param spacegroup the spacegroup.
554
\param cell the cell. */
555
MiniMol::MiniMol( const Spacegroup& spacegroup, const Cell& cell )
557
init( spacegroup, cell );
558
Message::message( message_ctor_mmodel );
561
/*! The object is initialised with no atoms.
562
\param spacegroup the spacegroup.
563
\param cell the cell. */
564
void MiniMol::init( const Spacegroup& spacegroup, const Cell& cell )
566
spacegroup_ = spacegroup;
570
bool MiniMol::is_null() const
571
{ return ( spacegroup_.is_null() || cell_.is_null() ); }
576
// UTILITY FUNCTIONS:
577
// e.g. protein specific tools.
580
/*! A carbonyl oxygen is added to this residue if the supplied residue
581
contains an appriate N atom bonded to the C. Otherwise, nothing
583
\param next The next monomer in the chain.
585
void MMonomer::protein_mainchain_build_carbonyl_oxygen( const MMonomer& next ) {
586
// check for mainchain atoms
587
int a1 = lookup( " CA ", MM::ANY );
588
int c1 = lookup( " C ", MM::ANY );
589
int n2 = next.lookup( " N ", MM::ANY );
590
if ( a1 < 0 || c1 < 0 || n2 < 0 ) return;
591
// get coordinates and check bonding
592
const clipper::Coord_orth cc1 = children[c1].coord_orth();
593
const clipper::Coord_orth ca1 = children[a1].coord_orth() - cc1;
594
const clipper::Coord_orth cn2 = next[n2].coord_orth() - cc1;
595
if ( cn2.lengthsq() > 2.2 ) return;
596
double uiso = children[c1].u_iso();
597
double occ = children[c1].occupancy();
598
// delete any existing O
599
int o1 = lookup( " O ", MM::ANY );
600
if ( o1 >= 0 ) children.erase( children.begin()+o1 );
602
const clipper::Vec3<> v0 = ca1.unit();
603
const clipper::Vec3<> v1 = ( clipper::Vec3<>::cross( v0, clipper::Vec3<>::cross( v0, cn2 ) ) ).unit();
604
MAtom atm = MAtom::null();
606
atm.set_element( "O" );
607
// length 1.24 angle 2.11
608
atm.set_coord_orth( cc1 + clipper::Coord_orth( -0.637*v0 + 1.064*v1 ) );
609
atm.set_occupancy( occ );
610
atm.set_u_iso( uiso );
614
// internal function for fast lookup in rotamer lib
615
int rotamer_find( String res, int rota ) {
616
if ( res.length() < 3 ) return 0;
618
int p2 = data::rotamer_data_size - 1;
619
while ( p2 - p1 > 1 ) {
620
int p = ( p1 + p2 ) / 2;
621
int s = strncmp( res.c_str(), data::rotamer_data[p].resname, 3 );
622
if ( s < 0 || ( s == 0 && rota <= data::rotamer_data[p].rota ) )
627
if ( strncmp( res.c_str(), data::rotamer_data[p2].resname, 3 ) == 0 &&
628
rota == data::rotamer_data[p2].rota ) return p2;
632
/*! \return The number of stored rotamers for this residue type.
634
int MMonomer::protein_sidechain_number_of_rotamers() const {
635
int r = rotamer_find( type(), 0 );
636
if ( r < 0 ) return 0;
637
return data::rotamer_data[r].num_rota;
641
\param n The number of the rotamer required.
642
\return The frequency of the given rotamer. */
643
ftype MMonomer::protein_sidechain_build_rotamer( const int& n ) {
644
int na = lookup( " CA ", MM::ANY );
645
int nc = lookup( " C ", MM::ANY );
646
int nn = lookup( " N ", MM::ANY );
647
if ( na < 0 || nc < 0 || nn < 0 ) return 0.0;
648
clipper::Coord_orth ca = children[na].coord_orth();
649
clipper::Coord_orth c1 = children[nc].coord_orth() - ca;
650
clipper::Coord_orth c2 = children[nn].coord_orth() - ca;
651
double uiso = children[na].u_iso();
652
double occ = children[na].occupancy();
653
// strip old side chain
654
for ( int i = children.size()-1; i >= 0; i-- )
655
if ( children[i].name() != " CA " && children[i].name() != " N " &&
656
children[i].name() != " C " && children[i].name() != " O " )
657
children.erase( children.begin()+i );
659
int r = rotamer_find( type(), n );
660
if ( r < 0 ) return 0.0;
661
if ( n >= data::rotamer_data[r].num_rota ) return -1.0;
662
// get rtop from standard orientation
663
const clipper::Vec3<> v1( (c1.unit()+c2.unit()).unit() );
664
const clipper::Vec3<> v2( clipper::Vec3<>::cross(c1,c2).unit() );
665
const clipper::Vec3<> v3( clipper::Vec3<>::cross(v1,v2).unit() );
666
const clipper::Mat33<> mr( v1[0], v2[0], v3[0],
668
v1[2], v2[2], v3[2] );
669
clipper::RTop_orth rtop( mr, ca );
671
MAtom atm = MAtom::null();
672
for ( int dr = 0; dr < data::rotamer_data[r].num_atom; dr++ ) {
674
String name = data::rotamer_data[i].atomname;
676
name = name.substr( 0, 2 );
678
atm.set_element( name );
679
atm.set_coord_orth( rtop * Coord_orth( data::rotamer_data[i].x,
680
data::rotamer_data[i].y,
681
data::rotamer_data[i].z ) );
682
atm.set_occupancy( occ );
683
atm.set_u_iso( uiso );
687
return data::rotamer_data[r].rota_prob;
691
/*! Test if the C of residue 1 is bonded to the N of residue 2,
692
within the distance r.
693
\param r1 The first residue.
694
\param r2 The second residue.
695
\param r The maximum allowed bond length.
696
\return true if N and C are present and bonded. */
697
bool MMonomer::protein_peptide_bond( const MMonomer& m1, const MMonomer& m2, ftype r ) {
698
int c1 = m1.lookup( " C ", MM::ANY );
699
int n2 = m2.lookup( " N ", MM::ANY );
700
if ( c1 >= 0 && n2 >= 0 )
701
if ( ( m1[c1].coord_orth() - m2[n2].coord_orth() ).lengthsq() < r*r )
706
/*! Return the Ramachadran angle in radians on -pi...pi.
707
To check the result, see clipper::Util::is_nan()
708
\param r1 The first residue.
709
\param r2 The second residue.
710
\return The torsion angle in radians, or NaN if atoms are missing. */
711
double MMonomer::protein_ramachandran_phi( const MMonomer& m1, const MMonomer& m2 )
713
ftype result = clipper::Util::nan();
714
int index_cx = m1.lookup( " C ", clipper::MM::ANY );
715
int index_n = m2.lookup( " N ", clipper::MM::ANY );
716
int index_ca = m2.lookup( " CA ", clipper::MM::ANY );
717
int index_c = m2.lookup( " C ", clipper::MM::ANY );
718
// if we have all three atoms, then add residue
719
if ( index_cx >= 0 && index_ca >= 0 && index_c >= 0 && index_n >= 0 ) {
720
Coord_orth coord_cx = m1[index_cx].coord_orth();
721
Coord_orth coord_n = m2[index_n ].coord_orth();
722
Coord_orth coord_ca = m2[index_ca].coord_orth();
723
Coord_orth coord_c = m2[index_c ].coord_orth();
725
result = Coord_orth::torsion( coord_cx, coord_n, coord_ca, coord_c );
730
/*! Return the Ramachadran angle in radians on -pi...pi.
731
To check the result, see clipper::Util::is_nan()
732
\param r1 The first residue.
733
\param r2 The second residue.
734
\return The torsion angle in radians, or NaN if atoms are missing. */
735
double MMonomer::protein_ramachandran_psi( const MMonomer& m1, const MMonomer& m2 )
737
ftype result = clipper::Util::nan();
738
int index_n = m1.lookup( " N ", clipper::MM::ANY );
739
int index_ca = m1.lookup( " CA ", clipper::MM::ANY );
740
int index_c = m1.lookup( " C ", clipper::MM::ANY );
741
int index_nx = m2.lookup( " N ", clipper::MM::ANY );
742
// if we have all three atoms, then add residue
743
if ( index_ca >= 0 && index_c >= 0 && index_n >= 0 && index_nx >= 0 ) {
744
Coord_orth coord_n = m1[index_n ].coord_orth();
745
Coord_orth coord_ca = m1[index_ca].coord_orth();
746
Coord_orth coord_c = m1[index_c ].coord_orth();
747
Coord_orth coord_nx = m2[index_nx].coord_orth();
749
result = Coord_orth::torsion( coord_n, coord_ca, coord_c, coord_nx );
755
} // namespace clipper
added
removed
clipper/minimol/minimol.cpp
