diff --git a/modelling/doc/index.rst b/modelling/doc/index.rst
index a0e0a7b513a7b97d8333033aa36de8c87fd85d0c..b7ae5aba1cc2dad4174bf0c68ea61ac2db61938d 100644
--- a/modelling/doc/index.rst
+++ b/modelling/doc/index.rst
@@ -106,9 +106,23 @@ Closing Gaps
.. autofunction:: CloseLargeDeletions
+.. function:: CountEnclosedGaps(mhandle, gap)
+ CountEnclosedInsertions(mhandle, gap)
+
+ Counts all gaps from `mhandle` which are fully enclosed by given `gap`.
+ This is either all gaps or only insertions.
+
+ :param mhandle: Modelling handle on which to apply change.
+ :type mhandle: :class:`ModellingHandle`
+ :param gap: Gap defining range in which gaps are to be removed.
+ :type gap: :class:`StructuralGap`
+
+ :return: Number of gaps.
+ :rtype: :class:`int`
+
.. function:: ClearGaps(mhandle, gap)
- Removes all gaps from mhandle which are fully enclosed by given gap.
+ Removes all gaps from `mhandle` which are fully enclosed by given `gap`.
:param mhandle: Modelling handle on which to apply change.
:type mhandle: :class:`ModellingHandle`
diff --git a/modelling/pymod/_closegaps.py b/modelling/pymod/_closegaps.py
index 8cbf086b15ad964e697de80e19b0997bcb2a6818..55f3ba293f9a08e0e72d8de247a9558865cb2589 100644
--- a/modelling/pymod/_closegaps.py
+++ b/modelling/pymod/_closegaps.py
@@ -502,6 +502,8 @@ def FillLoopsByDatabase(mhandle, scorer, fragment_db, structure_db,
# keep copy of original gap
gap_orig = mhandle.gaps[gap_idx].Copy()
+ if score_variant == 0:
+ gap_ins = CountEnclosedInsertions(mhandle, gap_orig)
# terminal gaps are no loops, database needs 2 stems to work on
if gap_orig.IsTerminal():
@@ -525,6 +527,9 @@ def FillLoopsByDatabase(mhandle, scorer, fragment_db, structure_db,
max_db_loop_len)
# iteratively extend actual_gap until we have enough loops
while actual_gap.length <= max_db_loop_len:
+ # check if we have enough
+ if found_loops >= max_loops_to_search:
+ break
# extend the gap
if not actual_extender.Extend():
break
@@ -532,10 +537,10 @@ def FillLoopsByDatabase(mhandle, scorer, fragment_db, structure_db,
if not actual_gap.before.IsValid() or \
not actual_gap.after.IsValid():
break
-
- # check if we have enough
- if found_loops >= max_loops_to_search:
- break
+ # skip if we try to merge insertions with score_variant=0
+ if score_variant == 0 and \
+ CountEnclosedInsertions(mhandle, actual_gap) != gap_ins:
+ continue
# get candidates for the current loop
candidates = loop.LoopCandidates.FillFromDatabase(
@@ -643,6 +648,8 @@ def FillLoopsByMonteCarlo(mhandle, scorer, torsion_sampler, max_loops_to_search=
# keep copy of original gap
gap_orig = mhandle.gaps[gap_idx].Copy()
+ if score_variant == 0:
+ gap_ins = CountEnclosedInsertions(mhandle, gap_orig)
# ignore terminal gaps
if gap_orig.IsTerminal():
@@ -678,6 +685,10 @@ def FillLoopsByMonteCarlo(mhandle, scorer, torsion_sampler, max_loops_to_search=
if not actual_gap.before.IsValid() or \
not actual_gap.after.IsValid():
break
+ # skip if we try to merge insertions with score_variant=0
+ if score_variant == 0 and \
+ CountEnclosedInsertions(mhandle, actual_gap) != gap_ins:
+ continue
# setup sampler, closer, scorer and cooler for MC
try:
diff --git a/modelling/pymod/export_model.cc b/modelling/pymod/export_model.cc
index 28e3b743f20410693dccfff250d6a8a1d6a793c3..5f457a2721aa48b1fa08dfd3e8ee52d86bb0090c 100644
--- a/modelling/pymod/export_model.cc
+++ b/modelling/pymod/export_model.cc
@@ -12,6 +12,14 @@ ModellingHandle (*BuildRawModelHandleList)(const ost::seq::AlignmentList&,
bool,
const String&)=&BuildRawModel;
+namespace {
+int WrapCountEnclosedGaps(ModellingHandle& mhandle, const StructuralGap& gap) {
+ return CountEnclosedGaps(mhandle, gap, false);
+}
+int WrapCountEnclosedIns(ModellingHandle& mhandle, const StructuralGap& gap) {
+ return CountEnclosedGaps(mhandle, gap, true);
+}
+} // anon ns
void export_model()
{
@@ -20,9 +28,11 @@ void export_model()
.def_readwrite("gaps", &ModellingHandle::gaps)
.def_readwrite("seqres", &ModellingHandle::seqres)
;
- def("ClearGaps",ClearGaps,(arg("gap")));
- def("MergeGaps",MergeGaps,(arg("index")));
- def("RemoveTerminalGaps",RemoveTerminalGaps);
+ def("ClearGaps", ClearGaps, (arg("mhandle"),arg("gap")));
+ def("CountEnclosedGaps", WrapCountEnclosedGaps, (arg("mhandle"),arg("gap")));
+ def("CountEnclosedInsertions", WrapCountEnclosedIns, (arg("mhandle"),arg("gap")));
+ def("MergeGaps", MergeGaps, (arg("mhandle"),arg("index")));
+ def("RemoveTerminalGaps", RemoveTerminalGaps, (arg("mhandle")));
def("BuildRawModel", BuildRawModelHandle,
(arg("aln"),
arg("include_ligands")=false,
diff --git a/modelling/src/model.cc b/modelling/src/model.cc
index c09654451e87e957fdee9f8a65fb5fec7bcfc298..4c99d1855da596fbb3f6007d8f30a80af6cd518d 100644
--- a/modelling/src/model.cc
+++ b/modelling/src/model.cc
@@ -48,6 +48,36 @@ bool CheckBackboneAtoms(ResidueView res)
}
+int CountEnclosedGaps(ModellingHandle& mhandle, const StructuralGap& gap,
+ bool insertions_only)
+{
+ int num_gaps = 0;
+
+ // extract info for this gap
+ String my_chain_name = gap.GetChainName();
+ ResNum my_b_res(0);
+ if (!gap.IsNTerminal()) my_b_res = gap.before.GetNumber();
+ ResNum my_a_res = my_b_res + gap.GetLength() + 1;
+
+ for (StructuralGapList::iterator i = mhandle.gaps.begin();
+ i != mhandle.gaps.end(); ++i) {
+
+ // skip different chains
+ if (i->GetChainName() != my_chain_name) continue;
+
+ // get res_nums
+ ResNum b_res(0);
+ if (!i->IsNTerminal()) b_res = i->before.GetNumber();
+ ResNum a_res = b_res + i->GetLength() + 1;
+
+ // check full overlap
+ if (b_res >= my_b_res && a_res <= my_a_res) {
+ if (!insertions_only || i->GetLength() > 0) ++num_gaps;
+ }
+ }
+ return num_gaps;
+}
+
int ClearGaps(ModellingHandle& mhandle, const StructuralGap& gap) {
// return -1 if no more gaps
diff --git a/modelling/src/model.hh b/modelling/src/model.hh
index 8e7be7d2ecc5b6fd0d5f3ed8894d65b4280d2273..6ee3516c4f87a49a8e49ac64dfd354e35fbc00a8 100644
--- a/modelling/src/model.hh
+++ b/modelling/src/model.hh
@@ -32,6 +32,10 @@ struct ModellingHandle {
ost::seq::SequenceList seqres;
};
+// see Python doc
+int CountEnclosedGaps(ModellingHandle& mhandle, const StructuralGap& gap,
+ bool insertions_only=false);
+
// see Python doc
int ClearGaps(ModellingHandle& mhandle, const StructuralGap& gap);
diff --git a/modelling/tests/test_modelling.py b/modelling/tests/test_modelling.py
index 909fc2db60254011ff2a181d154d22706132a827..cfe41f7c361bf06968c2ed1a24e890dd98077d35 100644
--- a/modelling/tests/test_modelling.py
+++ b/modelling/tests/test_modelling.py
@@ -129,6 +129,49 @@ class ModellingTests(unittest.TestCase):
with self.assertRaises(RuntimeError):
modelling.MergeGaps(mhandle, 0)
+ def testCountEnclosedGaps(self):
+ tpl = io.LoadPDB('data/2dbs.pdb')
+ aln = seq.CreateAlignment(
+ seq.CreateSequence('trg', 'GATLNGFTVPAGNTLVLN---PD--KG---ATVTMAGA'),
+ seq.CreateSequence('tpl', '--NGG--TL--LI--PNGTYHFLGIQMKSNVHIRVE--'))
+ aln.AttachView(1, tpl.CreateFullView())
+ mhandle = modelling.BuildRawModel(aln)
+ # check
+ self.assertEqual(len(mhandle.gaps), 8)
+ for g in mhandle.gaps:
+ self.assertEqual(modelling.CountEnclosedGaps(mhandle, g), 1)
+ self.assertEqual(modelling.CountEnclosedInsertions(mhandle, g),
+ 1 if (g.length > 0) else 0)
+
+ mych = mhandle.model.chains[0]
+ # none
+ mygap = modelling.StructuralGap(mych.FindResidue(15),
+ mych.FindResidue(17), "L")
+ self.assertEqual(modelling.CountEnclosedGaps(mhandle, mygap), 0)
+ self.assertEqual(modelling.CountEnclosedInsertions(mhandle, mygap), 0)
+ # extended singles
+ mygap = modelling.StructuralGap(mych.FindResidue(3),
+ mych.FindResidue(9), "LNGFT")
+ self.assertEqual(modelling.CountEnclosedGaps(mhandle, mygap), 1)
+ self.assertEqual(modelling.CountEnclosedInsertions(mhandle, mygap), 1)
+ mygap = modelling.StructuralGap(mych.FindResidue(20),
+ mych.FindResidue(22), "K")
+ self.assertEqual(modelling.CountEnclosedGaps(mhandle, mygap), 1)
+ self.assertEqual(modelling.CountEnclosedInsertions(mhandle, mygap), 0)
+ # doubles
+ mygap = modelling.StructuralGap(mych.FindResidue(3),
+ mych.FindResidue(12), "LNGFTVPA")
+ self.assertEqual(modelling.CountEnclosedGaps(mhandle, mygap), 2)
+ self.assertEqual(modelling.CountEnclosedInsertions(mhandle, mygap), 2)
+ mygap = modelling.StructuralGap(mych.FindResidue(13),
+ mych.FindResidue(19), "TLVLN")
+ self.assertEqual(modelling.CountEnclosedGaps(mhandle, mygap), 2)
+ self.assertEqual(modelling.CountEnclosedInsertions(mhandle, mygap), 1)
+ mygap = modelling.StructuralGap(mych.FindResidue(20),
+ mych.FindResidue(25), "KGAT")
+ self.assertEqual(modelling.CountEnclosedGaps(mhandle, mygap), 2)
+ self.assertEqual(modelling.CountEnclosedInsertions(mhandle, mygap), 0)
+
def testClearGapsExceptions(self):
tpl = io.LoadPDB('data/1mcg.pdb')
aln = seq.CreateAlignment(seq.CreateSequence('trg', 'DDFAGTHN'),