in layer2/Sculpt.cpp [390:1418]
void SculptMeasureObject(CSculpt * I, ObjectMolecule * obj, int state, int match_state,
int match_by_segment)
{
PyMOLGlobals *G = I->G;
int a, a0, a1, a2, a3, b0, b1, b2, b3, b4;
BondType *b;
float *v0, *v1, *v2, *v3, d, dummy;
CoordSet *cs;
int n0, n1, n2, n3;
int *planar = NULL;
int *linear = NULL;
int *single = NULL;
int *crdidx = NULL;
int nex = 1;
int *j, *k, xhash;
int ex_type;
AtomInfoType *ai, *ai1, *ai2, *obj_atomInfo;
int xoffset;
int use_cache = 1;
PRINTFD(G, FB_Sculpt)
" SculptMeasureObject-Debug: entered.\n" ENDFD;
if(match_state < 0)
match_state = state;
if(state < 0)
state = ObjectGetCurrentState(&obj->Obj, true);
ShakerReset(I->Shaker);
UtilZeroMem(I->NBHash, NB_HASH_SIZE * sizeof(int));
UtilZeroMem(I->EXHash, EX_HASH_SIZE * sizeof(int));
if((state >= 0) && (state < obj->NCSet) && (obj->CSet[state])) {
obj_atomInfo = obj->AtomInfo;
VLACheck(I->Don, int, obj->NAtom);
VLACheck(I->Acc, int, obj->NAtom);
ai = obj_atomInfo;
for(a = 0; a < obj->NAtom; a++) {
I->Don[a] = false;
I->Acc[a] = false;
AtomInfoCheckUniqueID(G, ai);
ai++;
}
ObjectMoleculeVerifyChemistry(obj, state);
ObjectMoleculeUpdateNeighbors(obj);
cs = obj->CSet[state];
use_cache = SettingGet_i(G, cs->Setting, obj->Obj.Setting, cSetting_sculpt_memory);
if(obj->NBond) {
int *neighbor = obj->Neighbor;
int n_atom = obj->NAtom;
planar = Alloc(int, n_atom);
linear = Alloc(int, n_atom);
single = Alloc(int, n_atom);
crdidx = Alloc(int, n_atom);
ai = obj_atomInfo;
for(a = 0; a < n_atom; a++) {
planar[a] = (ai->geom == cAtomInfoPlanar);
linear[a] = (ai->geom == cAtomInfoLinear);
single[a] = (ai->geom == cAtomInfoSingle);
if(obj->DiscreteFlag) {
if(cs == obj->DiscreteCSet[a]) {
a0 = obj->DiscreteAtmToIdx[a];
} else {
a0 = -1;
}
} else {
a0 = cs->AtmToIdx[a];
}
crdidx[a] = a0;
ai++;
}
/* brain-dead donor/acceptor assignment
* REPLACE later on with pattern-based system */
/* pass 1 */
b = obj->Bond;
for(a = 0; a < obj->NBond; a++) {
b1 = b->index[0];
b2 = b->index[1];
ai1 = obj_atomInfo + b1;
ai2 = obj_atomInfo + b2;
/* make blanket assumption that all nitrogens with
<3 bonds are donors -- we qualify this below... */
if(ai1->protons == cAN_N) {
n1 = neighbor[b1];
if(neighbor[n1] < 3) { /* N with L.P. */
I->Don[b1] = true;
}
}
if(ai2->protons == cAN_N) {
n2 = neighbor[b2];
if(neighbor[n2] < 3) { /* N with L.P. */
I->Don[b2] = true;
}
}
/* assume O is always an acceptor... */
if(ai1->protons == cAN_O)
I->Acc[b1] = true;
if(ai2->protons == cAN_O)
I->Acc[b2] = true;
b++;
}
/* pass 2 */
b = obj->Bond;
for(a = 0; a < obj->NBond; a++) {
b1 = b->index[0];
b2 = b->index[1];
/* nitrogens with lone pairs are acceptors
(not donors as assumed above) */
ai1 = obj_atomInfo + b1;
ai2 = obj_atomInfo + b2;
if(ai1->protons == cAN_N) {
if(b->order == 2) {
n1 = neighbor[b1];
if(neighbor[n1] < 3) { /* N with L.P. */
I->Acc[b1] = true;
I->Don[b1] = false;
}
}
}
if(ai2->protons == cAN_N) {
if(b->order == 2) {
n2 = neighbor[b2];
if(neighbor[n2] < 3) { /* N with L.P. */
I->Acc[b2] = true;
I->Don[b2] = false;
}
}
}
b++;
}
/* pass 3 */
b = obj->Bond;
for(a = 0; a < obj->NBond; a++) {
b1 = b->index[0];
b2 = b->index[1];
ai1 = obj_atomInfo + b1;
ai2 = obj_atomInfo + b2;
/* however, every NH is a donor,
even if it's SP2 */
if(ai1->protons == cAN_H) {
/* donors: any H attached to O, N */
switch (ai2->protons) {
case cAN_O:
I->Don[b1] = true;
I->Don[b2] = true; /* mark heavy atom too... */
break;
case cAN_N:
I->Don[b1] = true;
I->Don[b2] = true;
break;
}
} else if(ai2->protons == cAN_H) {
switch (ai1->protons) {
case cAN_O:
I->Don[b1] = true;
I->Don[b2] = true; /* mark heavy atom too... */
break;
case cAN_N:
I->Don[b1] = true;
I->Don[b2] = true; /* mark heavy atom too... */
break;
}
}
b++;
}
/* atom pass */
ai1 = obj_atomInfo;
for(a = 0; a < n_atom; a++) {
/* make sure all nonbonded atoms get categorized */
n0 = neighbor[a];
if(neighbor[n0] == 0) { /* nonbonded */
if(ai1->protons == cAN_O) {
I->Don[a] = true;
I->Acc[a] = true;
} else if(ai1->protons == cAN_N) {
I->Don[a] = true;
}
}
/*
if(I->Acc[a]) {
printf("ACC %s %s %s\n",ai1->chain,ai1->resi,ai1->name);
}
if(I->Don[a]) {
printf("DON %s %s %s\n",ai1->chain,ai1->resi,ai1->name);
} */
ai1++;
}
/* exclusions */
b = obj->Bond;
for(a = 0; a < obj->NBond; a++) {
b1 = b->index[0];
b2 = b->index[1];
ai1 = obj_atomInfo + b1;
ai2 = obj_atomInfo + b2;
xhash = ((b2 > b1) ? ex_hash(b1, b2) : ex_hash(b2, b1));
VLACheck(I->EXList, int, nex + 3);
j = I->EXList + nex;
*(j++) = *(I->EXHash + xhash);
if(b2 > b1) {
*(j++) = b1;
*(j++) = b2;
} else {
*(j++) = b2;
*(j++) = b1;
}
*(j++) = 2; /* 1-2 exclusion */
*(I->EXHash + xhash) = nex;
nex += 4;
a1 = crdidx[b1];
a2 = crdidx[b2];
if((a1 >= 0) && (a2 >= 0)) {
v1 = cs->Coord + 3 * a1;
v2 = cs->Coord + 3 * a2;
d = (float) diff3f(v1, v2);
if(use_cache) {
if(!SculptCacheQuery(G, cSculptBond,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b2].unique_id, 0, 0, &d))
SculptCacheStore(G, cSculptBond,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b2].unique_id, 0, 0, d);
}
ShakerAddDistCon(I->Shaker, b1, b2, d, cShakerDistBond, 1.0F);
/* NOTE: storing atom indices, not coord. ind.! */
}
b++;
}
/* triangle relationships */
{
ATLCall atl;
ai1 = obj_atomInfo;
atl.G = I->G;
atl.Shaker = I->Shaker;
atl.ai = obj_atomInfo;
atl.cSet = cs;
if(obj->DiscreteFlag) {
atl.atm2idx = obj->DiscreteAtmToIdx;
atl.discCSet = obj->DiscreteCSet;
} else {
atl.atm2idx = cs->AtmToIdx;
atl.discCSet = NULL;
}
atl.coord = cs->Coord;
atl.neighbor = neighbor;
atl.min = SettingGet_i(G, cs->Setting, obj->Obj.Setting, cSetting_sculpt_tri_min);
atl.max = SettingGet_i(G, cs->Setting, obj->Obj.Setting, cSetting_sculpt_tri_max);
atl.mode =
SettingGet_i(G, cs->Setting, obj->Obj.Setting, cSetting_sculpt_tri_mode);
for(a = 0; a < n_atom; a++) {
atl.atom0 = a;
/* clear the flag -- TODO replace with array */
{
int aa;
ai = obj_atomInfo;
for(aa = 0; aa < n_atom; aa++) {
ai->temp1 = false;
ai++;
}
}
ai1->temp1 = true;
add_triangle_limits(&atl, a, a, 0.0F, 1);
ai1++;
}
}
/* if we have a match state, establish minimum distances */
if((match_state >= 0) && (match_state < obj->NCSet) && (!obj->DiscreteFlag)) {
CoordSet *cs2 = obj->CSet[match_state];
int n_site = 0;
if(cs2) {
float minim_min =
SettingGet_f(G, cs->Setting, obj->Obj.Setting, cSetting_sculpt_min_min);
float minim_max =
SettingGet_f(G, cs->Setting, obj->Obj.Setting, cSetting_sculpt_min_max);
float maxim_min =
SettingGet_f(G, cs->Setting, obj->Obj.Setting, cSetting_sculpt_max_min);
float maxim_max =
SettingGet_f(G, cs->Setting, obj->Obj.Setting, cSetting_sculpt_max_max);
int *site = Calloc(int, n_atom);
float *weight = Calloc(float, n_atom);
/* first, find candidate atoms with sufficient connectivity */
CountCall cnt;
cnt.ai = obj_atomInfo;
cnt.neighbor = neighbor;
cnt.atm2idx1 = cs->AtmToIdx;
cnt.atm2idx2 = cs2->AtmToIdx;
{
int aa;
ai = obj_atomInfo;
for(aa = 0; aa < n_atom; aa++) {
ai->temp1 = false;
ai++;
}
}
ai1 = obj_atomInfo;
for(b0 = 0; b0 < n_atom; b0++) {
int n_qual_branch = 0, cb;
int adj_site = false;
ai1->temp1 = true;
n0 = neighbor[b0] + 1;
while((b1 = neighbor[n0]) >= 0) {
if(site[b1]) {
adj_site = true;
break;
}
cb = count_branch(&cnt, b1, 3);
if(cb > 3) {
n_qual_branch++;
}
n0 += 2;
}
ai1->temp1 = false;
if((n_qual_branch > 2) && (!adj_site)) {
site[b0] = 10;
} else if(!adj_site) {
const char * name = LexStr(G, ai1->name);
switch (name[0]) {
case 'O':
if(!name[1])
if(AtomInfoKnownPolymerResName(LexStr(G, ai1->resn)))
site[b0] = 40; /* main-chain carbonyl */
break;
case 'C':
switch (name[1]) {
case 'Z':
switch (name[2]) {
case 0:
if(ai1->resn == G->lex_const.ARG)
site[b0] = 20; /* ARG/CZ */
else if(ai1->resn == G->lex_const.TYR)
site[b0] = 20; /* TYR/CZ */
else if(ai1->resn == G->lex_const.PHE)
site[b0] = 20; /* PHE/CZ */
break;
}
break;
case 'E':
switch (name[2]) {
case 0:
if(ai1->resn == G->lex_const.LYS)
site[b0] = 20; /* LYS/CE */
break;
}
break;
case 'D':
switch (name[2]) {
case 0:
if(ai1->resn == G->lex_const.GLU)
site[b0] = 20; /* GLU/CD */
else if(ai1->resn == G->lex_const.GLN)
site[b0] = 20; /* GLN/CD */
break;
}
break;
case 'G':
switch (name[2]) {
case 0:
if(ai1->resn == G->lex_const.LEU)
site[b0] = 20; /* LEU/CG */
else if(ai1->resn == G->lex_const.ASP)
site[b0] = 20; /* ASP/CG */
else if(ai1->resn == G->lex_const.ASN)
site[b0] = 20; /* ASN/CG */
break;
}
break;
}
break;
case 'S':
switch (name[1]) {
case 'D':
switch (name[2]) {
case 0:
if(ai1->resn == G->lex_const.MET)
site[b0] = 20; /* MET/SD */
break;
}
break;
}
break;
}
}
ai1++;
}
for(b0 = 0; b0 < n_atom; b0++) {
if(site[b0]) {
weight[n_site] = 10.0F / site[b0];
site[n_site] = b0;
n_site++;
}
}
{
for(a0 = 0; a0 < n_site; a0++) {
for(a1 = a0 + 1; a1 < n_site; a1++) {
float wt = weight[a0] * weight[a1];
b0 = site[a0];
b1 = site[a1];
{
int i0a = cs->AtmToIdx[b0];
int i1a = cs->AtmToIdx[b1];
int i0b = cs2->AtmToIdx[b0];
int i1b = cs2->AtmToIdx[b1];
if((i0a >= 0) && (i1a >= 0) && (i0b >= 0) && (i1b >= 0) &&
((!match_by_segment)
|| (obj_atomInfo[b0].segi == obj_atomInfo[b1].segi))) {
float *v0a = cs->Coord + 3 * i0a;
float *v1a = cs->Coord + 3 * i1a;
float *v0b = cs2->Coord + 3 * i0b;
float *v1b = cs2->Coord + 3 * i1b;
float dist0, dist1, min_dist, max_dist;
dist0 = diff3f(v0a, v1a);
dist1 = diff3f(v0b, v1b);
min_dist = (dist0 < dist1) ? dist0 : dist1;
if((min_dist >= minim_min) && (min_dist <= minim_max)) {
ShakerAddDistCon(I->Shaker, b0, b1, min_dist, cShakerDistMinim, wt);
}
max_dist = (dist0 > dist1) ? dist0 : dist1;
if((max_dist >= maxim_min) && (max_dist <= maxim_max)) {
ShakerAddDistCon(I->Shaker, b0, b1, max_dist, cShakerDistMaxim, wt);
}
}
}
}
}
}
FreeP(weight);
FreeP(site);
}
}
/* now pick up those 1-3 interations */
/* b1-b0-b2 */
for(b0 = 0; b0 < n_atom; b0++) {
n0 = neighbor[b0] + 1;
while(neighbor[n0] >= 0) {
b1 = neighbor[n0];
n1 = n0 + 2;
while(neighbor[n1] >= 0) {
b2 = neighbor[n1];
xhash = ((b2 > b1) ? ex_hash(b1, b2) : ex_hash(b2, b1));
VLACheck(I->EXList, int, nex + 3);
j = I->EXList + nex;
*(j++) = *(I->EXHash + xhash);
if(b2 > b1) {
*(j++) = b1;
*(j++) = b2;
} else {
*(j++) = b2;
*(j++) = b1;
}
*(j++) = 3; /* 1-3 exclusion */
*(I->EXHash + xhash) = nex;
nex += 4;
a0 = crdidx[b0];
a1 = crdidx[b1];
a2 = crdidx[b2];
if((a0 >= 0) && (a1 >= 0) && (a2 >= 0)) {
v1 = cs->Coord + 3 * a1;
v2 = cs->Coord + 3 * a2;
d = (float) diff3f(v1, v2);
if(use_cache) {
if(!SculptCacheQuery(G, cSculptAngl,
obj_atomInfo[b0].unique_id,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b2].unique_id, 0, &d))
SculptCacheStore(G, cSculptAngl,
obj_atomInfo[b0].unique_id,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b2].unique_id, 0, d);
}
ShakerAddDistCon(I->Shaker, b1, b2, d, cShakerDistAngle, 1.0F);
if(linear[b0] && (linear[b1] || linear[b2])) {
if(use_cache) {
if(!SculptCacheQuery(G, cSculptLine,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b0].unique_id,
obj_atomInfo[b2].unique_id, 0, &dummy))
SculptCacheStore(G, cSculptLine,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b0].unique_id,
obj_atomInfo[b2].unique_id, 0, 0.0);
}
ShakerAddLineCon(I->Shaker, b1, b0, b2);
}
}
n1 += 2;
}
n0 += 2;
}
}
/* and record the pyramidal and planar geometries */
/* b1-b0-b2
* |
* b3 */
for(b0 = 0; b0 < n_atom; b0++) {
n0 = neighbor[b0] + 1;
while(neighbor[n0] >= 0) {
b1 = neighbor[n0];
n1 = n0 + 2;
while(neighbor[n1] >= 0) {
b2 = neighbor[n1];
n2 = n1 + 2;
while(neighbor[n2] >= 0) {
b3 = neighbor[n2];
a0 = crdidx[b0];
a1 = crdidx[b1];
a2 = crdidx[b2];
a3 = crdidx[b3];
if((a0 >= 0) && (a1 >= 0) && (a2 >= 0) && (a3 >= 0)) {
float d2 = 0.0F;
v0 = cs->Coord + 3 * a0;
v1 = cs->Coord + 3 * a1;
v2 = cs->Coord + 3 * a2;
v3 = cs->Coord + 3 * a3;
d = ShakerGetPyra(&d2, v0, v1, v2, v3);
if(fabs(d) < 0.05) {
planar[b0] = true;
}
if(planar[b0])
d = 0.0;
if(use_cache) {
if(!SculptCacheQuery(G, cSculptPyra,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b0].unique_id,
obj_atomInfo[b2].unique_id,
obj_atomInfo[b3].unique_id, &d))
SculptCacheStore(G, cSculptPyra,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b0].unique_id,
obj_atomInfo[b2].unique_id,
obj_atomInfo[b3].unique_id, d);
if(!SculptCacheQuery(G, cSculptPyra + 1,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b0].unique_id,
obj_atomInfo[b2].unique_id,
obj_atomInfo[b3].unique_id, &d2))
SculptCacheStore(G, cSculptPyra + 1,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b0].unique_id,
obj_atomInfo[b2].unique_id,
obj_atomInfo[b3].unique_id, d2);
}
ShakerAddPyraCon(I->Shaker, b0, b1, b2, b3, d, d2);
}
n2 += 2;
}
n1 += 2;
}
n0 += 2;
}
}
/* b1\b0_b2/b3 */
for(b0 = 0; b0 < n_atom; b0++) {
n0 = neighbor[b0] + 1;
while((b1 = neighbor[n0]) >= 0) {
n1 = neighbor[b0] + 1;
while((b2 = neighbor[n1]) >= 0) {
if(b1 != b2) {
n2 = neighbor[b2] + 1;
while((b3 = neighbor[n2]) >= 0) {
if((b3 != b0) && (b3 > b1)) {
if(!(planar[b0] || planar[b2] || linear[b0] || linear[b2])) {
int type;
if((obj_atomInfo[b0].protons == cAN_S) &&
(obj_atomInfo[b2].protons == cAN_S))
type = cShakerTorsDisulfide;
else
type = cShakerTorsSP3SP3;
ShakerAddTorsCon(I->Shaker, b1, b0, b2, b3, type);
}
if(planar[b0] && planar[b2]) {
/* special extra-rigid torsion for hydrogens on
planar acyclic systems (amides, etc.) */
if(((obj_atomInfo[b1].protons == cAN_H) && single[b1] &&
(obj_atomInfo[b3].protons != cAN_H) && planar[b3]) ||
((obj_atomInfo[b3].protons == cAN_H) && single[b3] &&
(obj_atomInfo[b1].protons != cAN_H) && planar[b1])) {
int cycle = 0;
/* b1\b0_b2/b3-b4-b5-b6-b7... */
int b5, b6, b7, b8, b9, b10;
int n4, n5, n6, n7, n8, n9;
n3 = neighbor[b2] + 1;
while((!cycle) && (b4 = neighbor[n3]) >= 0) {
if(b4 != b0) {
n4 = neighbor[b4] + 1;
while((!cycle) && (b5 = neighbor[n4]) >= 0) {
if(b5 != b2) {
n5 = neighbor[b5] + 1;
while((!cycle) && (b6 = neighbor[n5]) >= 0) {
if(b6 == b0) { /* 4-cycle */
cycle = 4;
} else if((b6 != b4) && (b6 != b2)) {
n6 = neighbor[b6] + 1;
while((!cycle) && (b7 = neighbor[n6]) >= 0) {
if(b7 == b0) { /* 5-cycle */
cycle = 5;
} else if((b7 != b5) && (b7 != b2)) {
n7 = neighbor[b7] + 1;
while((!cycle) && (b8 = neighbor[n7]) >= 0) {
if(b8 == b0) { /* 6-cycle */
cycle = 6;
} else if((b8 != b6) && (b8 != b2)) {
n8 = neighbor[b8] + 1;
while((!cycle) && (b9 = neighbor[n8]) >= 0) {
if(b9 == b0) { /* 7-cycle */
cycle = 7;
} else if((b9 != b7) && (b9 != b2)) {
n9 = neighbor[b9] + 1;
while((!cycle) && (b10 = neighbor[n9]) >= 0) {
if(b10 == b0) { /* 8-cycle */
cycle = 8;
}
n9 += 2;
}
}
n8 += 2;
}
}
n7 += 2;
}
}
n6 += 2;
}
}
n5 += 2;
}
}
n4 += 2;
}
}
n3 += 2;
}
if(!cycle) { /* don't add special amide constraints within small rings */
if(((obj_atomInfo[b1].protons == cAN_H) && single[b1] &&
(obj_atomInfo[b0].protons == cAN_N) &&
(obj_atomInfo[b2].protons == cAN_C) &&
(obj_atomInfo[b3].protons == cAN_O) && planar[b3]) ||
((obj_atomInfo[b1].protons == cAN_H) && single[b3] &&
(obj_atomInfo[b2].protons == cAN_N) &&
(obj_atomInfo[b0].protons == cAN_C) &&
(obj_atomInfo[b1].protons == cAN_O) && planar[b1])) {
/* biased, asymmetric term for amides */
ShakerAddTorsCon(I->Shaker, b1, b0, b2, b3, cShakerTorsAmide);
} else {
/* biased, symmetric term for all others */
ShakerAddTorsCon(I->Shaker, b1, b0, b2, b3, cShakerTorsFlat);
}
}
}
}
/* check 1-4 exclusion */
xhash = ex_hash(b1, b3);
ex_type = 4;
xoffset = *(I->EXHash + xhash);
while(xoffset) {
k = I->EXList + xoffset;
if((abs(*(k + 3)) == 4) && (*(k + 1) == b1) && (*(k + 2) == b3)) {
if((b0 != *(k + 4)) && (b2 != *(k + 5))) {
if(planar[b0] && planar[b2] &&
planar[*(k + 4)] && planar[*(k + 5)]) {
/* two planar paths -> likely a planar aromatic system */
*(k + 3) = -4;
}
}
ex_type = 0; /* duplicate, skip */
break;
}
xoffset = *k;
}
if(ex_type) {
VLACheck(I->EXList, int, nex + 5);
j = I->EXList + nex;
*(j++) = *(I->EXHash + xhash);
*(j++) = b1;
*(j++) = b3;
if(planar[b0] && planar[b2])
*(j++) = -4;
else
*(j++) = ex_type;
*(j++) = b0;
*(j++) = b2;
*(I->EXHash + xhash) = nex;
nex += 6;
}
/* planarity */
a0 = crdidx[b0];
a1 = crdidx[b1];
a2 = crdidx[b2];
a3 = crdidx[b3];
if((a0 >= 0) && (a1 >= 0) && (a2 >= 0) && (a3 >= 0)) {
v0 = cs->Coord + 3 * a0;
v1 = cs->Coord + 3 * a1;
v2 = cs->Coord + 3 * a2;
v3 = cs->Coord + 3 * a3;
d = 0.0;
if(planar[b0] && planar[b2]) {
float deg = get_dihedral3f(v1, v0, v2, v3);
if(fabs(deg) < deg_to_rad(10.0))
d = 1.0;
else if(fabs(deg) > deg_to_rad(170))
d = -1.0;
{
int cycle = false;
/* look for 4, 5, 6, 7, or 8 cycle that
connects back to b1 if found, then this
planar system is fixed (either at zero
or 180 -- it can't flip it over) */
/* b1\b0_b2/b3-b4-b5-b6-b7... */
int b5, b6, b7, b8, b9, b10;
int n4, n5, n6, n7, n8, n9;
n3 = neighbor[b2] + 1;
while((!cycle) && (b4 = neighbor[n3]) >= 0) {
if(b4 != b0) {
n4 = neighbor[b4] + 1;
while((!cycle) && (b5 = neighbor[n4]) >= 0) {
if(b5 != b2) {
n5 = neighbor[b5] + 1;
while((!cycle) && (b6 = neighbor[n5]) >= 0) {
if(b6 == b0) { /* 4-cycle */
cycle = 4;
} else if((b6 != b4) && (b6 != b2)) {
n6 = neighbor[b6] + 1;
while((!cycle) && (b7 = neighbor[n6]) >= 0) {
if(b7 == b0) { /* 5-cycle */
cycle = 5;
} else if((b7 != b5) && (b7 != b2)) {
n7 = neighbor[b7] + 1;
while((!cycle) && (b8 = neighbor[n7]) >= 0) {
if(b8 == b0) { /* 6-cycle */
cycle = 6;
} else if((b8 != b6) && (b8 != b2)) {
n8 = neighbor[b8] + 1;
while((!cycle) && (b9 = neighbor[n8]) >= 0) {
if(b9 == b0) { /* 7-cycle */
cycle = 7;
} else if((b9 != b7) && (b9 != b2)) {
n9 = neighbor[b9] + 1;
while((!cycle)
&& (b10 = neighbor[n9]) >= 0) {
if(b10 == b0) { /* 8-cycle */
cycle = 8;
}
n9 += 2;
}
}
n8 += 2;
}
}
n7 += 2;
}
}
n6 += 2;
}
}
n5 += 2;
}
}
n4 += 2;
}
}
n3 += 2;
}
/* don't get jacked by pseudo-planar PRO */
if(((obj_atomInfo[b0].protons != cAN_N) ||
(!WordMatchExact(G, obj_atomInfo[b0].resn, G->lex_const.PRO, true))) &&
((obj_atomInfo[b2].protons != cAN_N) ||
(!WordMatchExact(G, obj_atomInfo[b2].resn, G->lex_const.PRO, true)))) {
if(use_cache) {
if(!SculptCacheQuery(G, cSculptPlan,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b0].unique_id,
obj_atomInfo[b2].unique_id,
obj_atomInfo[b3].unique_id, &d))
SculptCacheStore(G, cSculptPlan,
obj_atomInfo[b1].unique_id,
obj_atomInfo[b0].unique_id,
obj_atomInfo[b2].unique_id,
obj_atomInfo[b3].unique_id, d);
}
ShakerAddPlanCon(I->Shaker, b1, b0, b2, b3, d, cycle);
if(planar[b1] && planar[b3] && ((cycle == 5) || (cycle == 6))) {
/* also add minimum distance constraints to keep small rings from folding */
d = (float) diff3f(v1, v3);
ShakerAddDistCon(I->Shaker, b1, b3, d, cShakerDistBond, 1.0F);
}
}
}
}
}
}
n2 += 2;
}
}
n1 += 2;
}
n0 += 2;
}
}
/* add 1,5 exclusions for hydrogens off arg-like planar systems */
/* b1\b0_b2_b3/b4 */
for(b0 = 0; b0 < n_atom; b0++) {
n0 = neighbor[b0] + 1;
while((b1 = neighbor[n0]) >= 0) {
if(obj_atomInfo[b1].protons == cAN_H) {
n1 = neighbor[b0] + 1;
while((b2 = neighbor[n1]) >= 0) {
if(b1 != b2) {
n2 = neighbor[b2] + 1;
while((b3 = neighbor[n2]) >= 0) {
if(b3 != b0) {
if(planar[b0] && planar[b2] && planar[b3]) {
n3 = neighbor[b3] + 1;
while((b4 = neighbor[n3]) >= 0) {
if((b4 != b2) && (b4 > b1) && (obj_atomInfo[b4].protons == cAN_H)) {
xhash = ex_hash(b1, b4);
ex_type = 5;
xoffset = *(I->EXHash + xhash);
while(xoffset) {
k = I->EXList + xoffset;
if(((*(k + 3)) == ex_type) &&
(*(k + 1) == b1) && (*(k + 2) == b4)) {
ex_type = 0; /* duplicate, skip */
break;
}
xoffset = *k;
}
if(ex_type) {
VLACheck(I->EXList, int, nex + 6);
j = I->EXList + nex;
*(j++) = *(I->EXHash + xhash);
*(j++) = b1;
*(j++) = b4;
*(j++) = ex_type;
*(j++) = b0;
*(j++) = b2;
*(j++) = b3;
*(I->EXHash + xhash) = nex;
nex += 6;
}
}
n3 += 2;
}
}
}
n2 += 2;
}
}
n1 += 2;
}
}
n0 += 2;
}
}
{
/* longer-range exclusions (1-5,1-6,1-7,1-8,1-9) -- only locate & store when needed */
int mask =
SettingGet_i(G, cs->Setting, obj->Obj.Setting, cSetting_sculpt_field_mask);
int max_excl =
SettingGet_i(G, cs->Setting, obj->Obj.Setting, cSetting_sculpt_avd_excl);
if(max_excl > 9)
max_excl = 9;
if((cSculptAvoid & mask) && (max_excl > 4)) {
int b_stack[10];
int n_stack[10];
int stop_depth = max_excl - 1;
int depth;
int bd, skip;
for(b0 = 0; b0 < n_atom; b0++) {
b_stack[0] = b0;
n_stack[0] = neighbor[b_stack[0]] + 1;
depth = 0;
while(depth >= 0) {
if((bd = neighbor[n_stack[depth]]) < 0) {
depth--;
if(depth >= 0) { /* iterate next atom */
n_stack[depth] += 2;
}
} else {
skip = (depth == stop_depth);
if(!skip) {
for(a = 0; a < depth; a++) {
if(b_stack[a] == bd) {
skip = true;
break;
}
}
}
if(!skip) {
depth++;
b_stack[depth] = bd;
n_stack[depth] = neighbor[bd] + 1;
if((depth > 3) && (b0 < bd)) {
xhash = ex_hash(b0, bd);
VLACheck(I->EXList, int, nex + 3);
j = I->EXList + nex;
*(j++) = *(I->EXHash + xhash);
*(j++) = b0;
*(j++) = bd;
*(j++) = depth + 1; /* 1-5, 1-6, 1-7 etc. */
*(I->EXHash + xhash) = nex;
nex += 4;
}
} else {
n_stack[depth] += 2;
}
}
}
}
}
}
FreeP(planar);
FreeP(linear);
FreeP(single);
FreeP(crdidx);
}
}
PRINTFB(G, FB_Sculpt, FB_Blather)
" Sculpt: I->Shaker->NDistCon %d\n", I->Shaker->NDistCon ENDFB(G);
PRINTFB(G, FB_Sculpt, FB_Blather)
" Sculpt: I->Shaker->NPyraCon %d\n", I->Shaker->NPyraCon ENDFB(G);
PRINTFB(G, FB_Sculpt, FB_Blather)
" Sculpt: I->Shaker->NPlanCon %d\n", I->Shaker->NPlanCon ENDFB(G);
PRINTFD(G, FB_Sculpt)
" SculptMeasureObject-Debug: leaving...\n" ENDFD;
}