From 6ae59864a6f04e79bc98fdb445ef0a50f4c1a39a Mon Sep 17 00:00:00 2001 From: Manuel Lera-Ramirez Date: Thu, 7 Nov 2024 11:11:52 +0000 Subject: [PATCH] closes #279 again (for circular seqs) (#324) --- src/pydna/amplify.py | 6 ++++-- tests/test_module_amplify.py | 5 +++++ 2 files changed, 9 insertions(+), 2 deletions(-) diff --git a/src/pydna/amplify.py b/src/pydna/amplify.py index 3c12854a..09ff68e5 100644 --- a/src/pydna/amplify.py +++ b/src/pydna/amplify.py @@ -351,11 +351,13 @@ def products(self): feats = self.template[ fp.position - fp._fp : rp.position + rp._fp ].features # Save features covered by primers - shift_amount = len(fp.tail) - feats = [_shift_feature(f, shift_amount, None) for f in feats] tpl = self.template else: continue + # Shift features to the right if there was a tail + shift_amount = len(fp.tail) + feats = [_shift_feature(f, shift_amount, None) for f in feats] + if tpl.circular and fp.position == rp.position: prd = _Dseqrecord(fp) + _Dseqrecord(rp).reverse_complement() else: diff --git a/tests/test_module_amplify.py b/tests/test_module_amplify.py index eeb4d895..8aa869d0 100644 --- a/tests/test_module_amplify.py +++ b/tests/test_module_amplify.py @@ -799,6 +799,11 @@ def test_annotation(): assert pcr_product.features[0].location.extract(pcr_product).seq == dsr.seq + # Also works in circular sequences + dsr_circ = dsr.looped() + pcr_product_circ = pcr(forward_primer, reverse_primer, dsr_circ) + assert str(pcr_product_circ.features[0].location.extract(pcr_product_circ).seq) == str(dsr_circ.seq) + if __name__ == "__main__": pytest.main([__file__, "-vv", "-s"])