</div>
</div>
+<div id="method-M000261" class="method-detail">
+<a name="find_close_atoms"></a>
+<div class="method-heading">
+<span class="method-name">find_close_atoms(atom, limit = 1.2) → Array of Integers<br />
+</span>
+</div>
+<div class="method-description">
+<p>
+Find atoms that are within the threshold distance from the given atom. If limit is a positive number, the threshold distance is the sum of the vdw radii times limit. If limit is a negative number, its absolute value is used for the threshold distance in angstrom. If limit is not given, a default value of 1.2 is used. The number of the newly created bonds is returned. An array of atom indices is returned. If no atoms are found, an empty array is returned.
+</p>
+</div>
+</div>
+
<div id="method-M000294" class="method-detail">
<a name="M000294"></a>
<div class="method-heading">
</div>
<div class="method-description">
<p>
-Create bonds between atoms that are close enough,
-i.e. the interatomic distance is smaller than the sum of the vdw radii
-times the argument limit. Returns the number of the newly created bonds.
-This operation is undoable.
+Create bonds between atoms that are within the threshold distance. If limit is a positive number, the threshold distance is the sum of the vdw radii times limit. If limit is a negative number, its absolute value is used for the threshold distance in angstrom. If limit is not given, a default value of 1.2 is used. The number of the newly created bonds is returned. This operation is undoable.
</p>
</div>
</div>
return sMoleculeSearchAcrossAtomGroup(mp->nimpropers, mp->impropers, 4, atomgroup, "impropers");
}
+/* Subroutine for MoleculeGuessBonds. It can be also used independently, but make sure that *outNbonds/*outBonds
+ _correctly_ represents an array of two integers (as in mp->nbonds/mp->bonds). */
+/* Find atoms within the given "distance" from the given atom. */
+/* If limit is negative, its absolute value denotes the threshold distance in angstrom; otherwise,
+ the threshold distance is given by the sum of van der Waals radii times limit. */
+/* If triangle is non-zero, then only atoms with lower indexes than index are looked for. */
+int
+MoleculeFindCloseAtoms(Molecule *mp, Int index, Double limit, Int *outNbonds, Int **outBonds, Int triangle)
+{
+ Int n1, n2, j, nlim, newbond[2];
+ Double a1, a2, alim;
+ Vector dr, r1, r2;
+ Atom *ap = ATOM_AT_INDEX(mp->atoms, index);
+ n1 = ap->atomicNumber;
+ if (n1 >= 0 && n1 < gCountElementParameters)
+ a1 = gElementParameters[n1].radius;
+ else a1 = gElementParameters[6].radius;
+ r1 = ap->r;
+ nlim = (triangle ? index : mp->natoms);
+ for (j = 0; j < nlim; j++) {
+ Atom *bp = ATOM_AT_INDEX(mp->atoms, j);
+ if (index == j)
+ continue;
+ n2 = bp->atomicNumber;
+ if (n2 >= 0 && n2 < gCountElementParameters)
+ a2 = gElementParameters[n2].radius;
+ else a2 = gElementParameters[6].radius;
+ r2 = bp->r;
+ VecSub(dr, r1, r2);
+ if (limit < 0)
+ alim = -limit;
+ else
+ alim = limit * (a1 + a2);
+ if (VecLength2(dr) < alim * alim) {
+ newbond[0] = index;
+ newbond[1] = j;
+ /* MoleculeAddBonds(mp, 1, newbonds); */
+ AssignArray(outBonds, outNbonds, sizeof(Int) * 2, *outNbonds, newbond);
+ }
+ }
+ return 0;
+}
+
/* Guess the bonds from the coordinates */
+/* If limit is negative, its absolute value denotes the threshold distance in angstrom; otherwise,
+ the threshold distance is given by the sum of van der Waals radii times limit. */
int
MoleculeGuessBonds(Molecule *mp, Double limit, Int *outNbonds, Int **outBonds)
{
- int i, j, n1, n2;
- Int nbonds, *bonds;
+ Int nbonds, *bonds, i, newbond[2];
+/* int i, j, n1, n2;
Atom *ap, *bp;
Vector r1, r2, dr;
Double a1, a2, alim;
Int newbond[2];
- ElementPar *p = gElementParameters;
+ ElementPar *p = gElementParameters; */
nbonds = 0;
bonds = NULL;
for (i = 0; i < mp->natoms; i++) {
+ MoleculeFindCloseAtoms(mp, i, limit, &nbonds, &bonds, 1);
+ /*
ap = ATOM_AT_INDEX(mp->atoms, i);
n1 = ap->atomicNumber;
if (n1 >= 0 && n1 < gCountElementParameters)
a1 = p[n1].radius;
else a1 = p[6].radius;
r1 = ap->r;
- /* if (mp->is_xtal_coord)
- TransformVec(&r1, mp->cell->tr, &r1); */
for (j = 0; j < i; j++) {
bp = ATOM_AT_INDEX(mp->atoms, j);
n2 = bp->atomicNumber;
a2 = p[n2].radius;
else a2 = p[6].radius;
r2 = bp->r;
- /* if (mp->is_xtal_coord)
- TransformVec(&r2, mp->cell->tr, &r2); */
VecSub(dr, r1, r2);
- alim = limit * (a1 + a2);
+ if (limit < 0)
+ alim = -limit;
+ else
+ alim = limit * (a1 + a2);
if (VecLength2(dr) < alim * alim) {
newbond[0] = i;
newbond[1] = j;
- /* MoleculeAddBonds(mp, 1, newbonds); */
AssignArray(&bonds, &nbonds, sizeof(Int) * 2, nbonds, newbond);
}
}
+ */
}
if (nbonds > 0) {
newbond[0] = kInvalidIndex;
int MoleculeAnalyzeAtomName(const char *s, char *resName, int *resSeq, char *atomName);
int MoleculeAtomIndexFromString(Molecule *mp, const char *s);
+int MoleculeFindCloseAtoms(Molecule *mp, Int index, Double limit, Int *outNbonds, Int **outBonds, Int triangle);
int MoleculeGuessBonds(Molecule *mp, Double limit, Int *outNbonds, Int **outBonds);
int MoleculeAreAtomsConnected(Molecule *mp, int n1, int n2);
int MoleculeRebuildTablesFromConnects(Molecule *mp);
/*
* call-seq:
+ * find_close_atoms(atom, limit = 1.2) -> array of Integers (atom indices)
+ *
+ * Find atoms that are within the threshold distance from the given atom.
+ * If limit is a positive number, the threshold distance is the sum of the vdw radii times limit.
+ * If limit is a negative number, its absolute value is used for the threshold distance in angstrom.
+ * If limit is not given, a default value of 1.2 is used.
+ * An array of atom indices is returned. If no atoms are found, an empty array is returned.
+ */
+static VALUE
+s_Molecule_FindCloseAtoms(int argc, VALUE *argv, VALUE self)
+{
+ Molecule *mol;
+ VALUE aval, limval;
+ double limit;
+ Int n1, nbonds, *bonds;
+ Data_Get_Struct(self, Molecule, mol);
+ rb_scan_args(argc, argv, "11", &aval, &limval);
+ n1 = s_Molecule_AtomIndexFromValue(mol, aval);
+ if (limval == Qnil)
+ limit = 1.2;
+ else
+ limit = NUM2DBL(rb_Float(limval));
+ nbonds = 0; /* This initialization is necessary: see comments in MoleculeFindCloseAtoms() */
+ bonds = NULL;
+ MoleculeFindCloseAtoms(mol, n1, limit, &nbonds, &bonds, 0);
+ aval = rb_ary_new();
+ if (nbonds > 0) {
+ for (n1 = 0; n1 < nbonds; n1++)
+ rb_ary_push(aval, INT2NUM(bonds[n1 * 2 + 1]));
+ free(bonds);
+ }
+ return aval;
+}
+
+/*
+ * call-seq:
* guess_bonds(limit = 1.2) -> Integer
*
- * Create bonds between atoms that are 'close enough', i.e. the interatomic distance is
- * smaller than the sum of the vdw radii times the argument 'limit'. If limit is not
- * given, a default value of 1.2 is used.
+ * Create bonds between atoms that are within the threshold distance.
+ * If limit is a positive number, the threshold distance is the sum of the vdw radii times limit.
+ * If limit is a negative number, its absolute value is used for the threshold distance in angstrom.
+ * If limit is not given, a default value of 1.2 is used.
* The number of the newly created bonds is returned.
* This operation is undoable.
*/
rb_define_method(rb_cMolecule, "assign_residue", s_Molecule_AssignResidue, 2);
rb_define_method(rb_cMolecule, "offset_residue", s_Molecule_OffsetResidue, 2);
rb_define_method(rb_cMolecule, "renumber_atoms", s_Molecule_RenumberAtoms, 1);
+ rb_define_method(rb_cMolecule, "find_close_atoms", s_Molecule_FindCloseAtoms, -1);
rb_define_method(rb_cMolecule, "guess_bonds", s_Molecule_GuessBonds, -1);
rb_define_method(rb_cMolecule, "selection", s_Molecule_Selection, 0);
rb_define_method(rb_cMolecule, "selection=", s_Molecule_SetSelection, 1);
OSDN Git Service
