diff --git a/modules/mol/base/pymod/export_entity_view.cc b/modules/mol/base/pymod/export_entity_view.cc
index 5e3dc12d2a9c15f10d3c4b9ce0df9a922554d16b..55d408365dc1c982c5370034af96a7f2f98de78b 100644
--- a/modules/mol/base/pymod/export_entity_view.cc
+++ b/modules/mol/base/pymod/export_entity_view.cc
@@ -20,6 +20,8 @@
#include <boost/python/suite/indexing/vector_indexing_suite.hpp>
using namespace boost::python;
+#include <ost/export_helper/pair_to_tuple_conv.hh>
+
#include <ost/mol/entity_view.hh>
#include <ost/mol/query.hh>
#include <ost/mol/mol.hh>
@@ -184,6 +186,8 @@ void export_EntityView()
.def("GetBounds", &EntityView::GetBounds)
.add_property("bounds", &EntityView::GetBounds)
;
+ to_python_converter<std::pair<EntityView, EntityView>,
+ PairToTupleConverter<EntityView, EntityView> >();
def("Union", &Union);
def("Difference", &Difference);
def("Intersection", &Intersection);
diff --git a/modules/mol/base/src/entity_view.hh b/modules/mol/base/src/entity_view.hh
index 5d2a0848ced01cde8b11a0d5ad771d75d7c21008..6e5b2371695064033bfa4637d7a5e533e6c9957b 100644
--- a/modules/mol/base/src/entity_view.hh
+++ b/modules/mol/base/src/entity_view.hh
@@ -335,6 +335,8 @@ private:
EntityViewDataPtr data_;
};
+typedef std::pair<EntityView, EntityView> EntityViewPair;
+
}} // ns
#endif // OST_ENTITY_VIEW_HH
diff --git a/modules/seq/base/doc/seq.rst b/modules/seq/base/doc/seq.rst
index ade28f583e7b10f8dd9d7e0a4337e52347480b0b..f42a2f574c7afff92265932553bbe3c6d85f96b9 100644
--- a/modules/seq/base/doc/seq.rst
+++ b/modules/seq/base/doc/seq.rst
@@ -357,3 +357,22 @@ an alignment:
If master is set to -1, all sequences in the region are affected, otherwise
only the sequence at index equal to master is shifted.
+
+ .. method:: GetMatchingBackboneViews(index1=0, index2=1)
+
+ Returns a tuple of entity views containing matching backbone atoms for the
+ two sequences at index1 and index2, respectively. For each aligned column in
+ the alignment, backbone atoms are added to the view if both aligned residues
+ have them. It is guaranteed that the two views contain the same number of
+ atoms and that the order of the atoms in the two views is the same.
+
+ The output of this function can be used to superpose two structures with
+ :func:`~ost.mol.alg.SuperposeSVD`.
+
+
+ :param index1: The index of the first sequence
+
+ :param index2: The index of the second sequence.
+
+ :raises: In case one of the two sequences doesn't have an attached view, a
+ :exc:`RuntimeError` is raised.
diff --git a/modules/seq/base/pymod/const_seq_list_export_def.hh b/modules/seq/base/pymod/const_seq_list_export_def.hh
index 159857c3b9c0aed85b93783f07f0f13aaef7f101..cd1f52853b4ad241301d4c55b7e08752027a18cd 100644
--- a/modules/seq/base/pymod/const_seq_list_export_def.hh
+++ b/modules/seq/base/pymod/const_seq_list_export_def.hh
@@ -26,6 +26,7 @@
.def("IsValid", &C::IsValid) \
.def("Take", &C::Take) \
.def("Slice", &C::Slice) \
+ .def("FindSequence", &C::FindSequence) \
.def("SequencesHaveEqualLength", \
&C::SequencesHaveEqualLength) \
.def("__getitem__", &C::operator[]) \
diff --git a/modules/seq/base/pymod/export_sequence.cc b/modules/seq/base/pymod/export_sequence.cc
index 454245778ff6b6ca0b03dfa394d4830cb50e3394..c431dfa9b51a47a3c7a6a43eae25ab8a8c574a5e 100644
--- a/modules/seq/base/pymod/export_sequence.cc
+++ b/modules/seq/base/pymod/export_sequence.cc
@@ -21,8 +21,6 @@
#include <boost/python/register_ptr_to_python.hpp>
#include <boost/python/suite/indexing/vector_indexing_suite.hpp>
-
-#include <ost/export_helper/pair_to_tuple_conv.hh>
#include <ost/generic_property.hh>
#include <ost/export_helper/generic_property_def.hh>
#include <ost/info/info.hh>
@@ -299,6 +297,8 @@ void export_sequence()
.def("Copy", &AlignmentHandle::Copy)
.def("ToString", &AlignmentHandle::ToString)
.def("GetLength", &AlignmentHandle::GetLength)
+ .def("GetMatchingBackboneViews", &AlignmentHandle::GetMatchingBackboneViews,
+ (arg("index1")=0, arg("index2")=1))
.def("__len__", &AlignmentHandle::GetLength)
.def("GetSequences", &AlignmentHandle::GetSequences)
.def("GetCoverage", &AlignmentHandle::GetSequences)
@@ -357,8 +357,7 @@ void export_sequence()
.def("__getitem__", &do_slice_b)
;
implicitly_convertible<SequenceList, ConstSequenceList>();
- to_python_converter<std::pair<mol::EntityView, mol::EntityView>,
- PairToTupleConverter<mol::EntityView, mol::EntityView> >();
+
def("CreateSequenceList", &CreateSequenceList);
def("SequenceFromChain", seq_from_chain_a);
def("SequenceFromChain", seq_from_chain_b);
diff --git a/modules/seq/base/src/alignment_handle.cc b/modules/seq/base/src/alignment_handle.cc
index c39856e939fe0245c4c7c3ca33fe56af90ad5b8c..3e4ed4256f05bcef1128cdccef5585a3b1d8e582 100644
--- a/modules/seq/base/src/alignment_handle.cc
+++ b/modules/seq/base/src/alignment_handle.cc
@@ -20,8 +20,12 @@
/*
Author: Marco Biasini
*/
+
#include <ost/invalid_handle.hh>
+#include <ost/mol/residue_view.hh>
+#include <ost/mol/atom_view.hh>
#include <ost/seq/alignment_handle.hh>
+#include <ost/mol/residue_view.hh>
#include <ost/seq/impl/sequence_list_impl.hh>
#include <ost/seq/impl/sequence_impl.hh>
#include <ost/seq/sequence_list.hh>
@@ -270,5 +274,38 @@ Real AlignmentHandle::GetCoverage(int seq_index) const
return impl_->GetCoverage(seq_index);
}
+mol::EntityViewPair AlignmentHandle::GetMatchingBackboneViews(int idx0, int idx1) const
+{
+ this->CheckValidity();
+ const impl::SequenceImpl& s1=*impl_->GetSequence(idx0).get();
+ const impl::SequenceImpl& s2=*impl_->GetSequence(idx1).get();
+ if (!s1.HasAttachedView() || !s2.HasAttachedView()) {
+ throw std::runtime_error("both sequences must have a view attached");
+ }
+ mol::EntityView v1=s1.GetAttachedView().CreateEmptyView();
+ mol::EntityView v2=s2.GetAttachedView().CreateEmptyView();
+ for (int i=0; i<s1.GetLength(); ++i) {
+ if (s1[i]=='-' && s2[i]=='-') {
+ continue;
+ }
+ mol::ResidueView r1=s1.GetResidue(i);
+ mol::ResidueView r2=s2.GetResidue(i);
+ if (!r1.IsValid() || !r2.IsValid()) {
+ continue;
+ }
+ const char* bb_anames[]={"N", "CA", "C", "O"};
+ //for (size_t j=0; )
+ for (size_t j=0; j<4; ++j) {
+ mol::AtomView a1=r1.FindAtom(bb_anames[j]);
+ mol::AtomView a2=r2.FindAtom(bb_anames[j]);
+ if (!a1.IsValid() || !a2.IsValid()) {
+ continue;
+ }
+ v1.AddAtom(a1);
+ v2.AddAtom(a2);
+ }
+ }
+ return mol::EntityViewPair(v1, v2);
+}
}}
diff --git a/modules/seq/base/src/alignment_handle.hh b/modules/seq/base/src/alignment_handle.hh
index 151506d25042f87b56e022a40076b72ad11e7822..9a109215b7236dd7ccd3b7a9acf9ea838df48532 100644
--- a/modules/seq/base/src/alignment_handle.hh
+++ b/modules/seq/base/src/alignment_handle.hh
@@ -107,6 +107,8 @@ public:
void AttachView(int seq_index, const mol::EntityView& view,
const String& chain_name);
+
+ mol::EntityViewPair GetMatchingBackboneViews(int idx0=0, int idx1=1) const;
/// \brief set name of sequence
void SetSequenceName(int seq_index, const String& name);
diff --git a/modules/seq/base/src/impl/sequence_impl.hh b/modules/seq/base/src/impl/sequence_impl.hh
index fd7fd0dfe0bdea729a8de1402c7ed0b78ae39dde..beae61e03033f4820f36b2c062b6368cf155dd0e 100644
--- a/modules/seq/base/src/impl/sequence_impl.hh
+++ b/modules/seq/base/src/impl/sequence_impl.hh
@@ -137,6 +137,15 @@ public:
void Append(char olc);
+ char& operator[](size_t index)
+ {
+ return seq_string_[index];
+ }
+ char operator[](size_t index) const
+ {
+ return seq_string_[index];
+ }
+
private:
/// \brief Recalculates gap shifts from sequence.
@@ -160,7 +169,7 @@ private:
std::list<Shift> shifts_;
bool editing_;
int offset_;
- mol::EntityView attached_view_;
+ mol::EntityView attached_view_;
};
/// \internal
diff --git a/modules/seq/base/src/sequence_list.cc b/modules/seq/base/src/sequence_list.cc
index 699a7a0eb1d8c2e80d1f6b36aba790da3be44292..eb45568fdbb639c1a14eab4732f076f9d500bcee 100644
--- a/modules/seq/base/src/sequence_list.cc
+++ b/modules/seq/base/src/sequence_list.cc
@@ -72,6 +72,12 @@ SequenceList::Iterator SequenceList::End() const
return SequenceList::Iterator(impl_->End(), impl_->End());
}
+SequenceHandle SequenceList::FindSequence(const String& name) const
+{
+ this->CheckValidity();
+ return SequenceHandle(impl_->FindSequence(name));
+}
+
bool SequenceList::IsValid() const
{
return impl_.get()!=NULL;
@@ -226,6 +232,12 @@ ConstSequenceList ConstSequenceList::Slice(int first, int n) const
return ConstSequenceList(impl_->Slice(first, n));
}
+ConstSequenceHandle ConstSequenceList::FindSequence(const String& name) const
+{
+ this->CheckValidity();
+ return ConstSequenceHandle(impl_->FindSequence(name));
+}
+
ConstSequenceList DLLEXPORT_OST_SEQ CreateConstSequenceList()
{
return ConstSequenceList(impl::SequenceListImplPtr(new impl::SequenceListImpl));
diff --git a/modules/seq/base/src/sequence_list.hh b/modules/seq/base/src/sequence_list.hh
index 62eea1041311f95fc697144725948883a211f9c7..71067438d0404a0af7c41a6138b3fb09c403433f 100644
--- a/modules/seq/base/src/sequence_list.hh
+++ b/modules/seq/base/src/sequence_list.hh
@@ -61,6 +61,8 @@ public:
/// \brief create a sequence list from the given splice interval
ConstSequenceList Slice(int first, int n) const;
+
+ ConstSequenceHandle FindSequence(const String& name) const;
int GetMinLength() const;
int GetMaxLength() const;
/// \internal
@@ -107,6 +109,7 @@ public:
int GetMinLength() const;
int GetMaxLength() const;
+ SequenceHandle FindSequence(const String& name) const;
// \internal
impl::SequenceListImplPtr& Impl() const;
SequenceList(const impl::SequenceListImplPtr& impl);