From 465d369a9c8c58b0a291b6582a884d3c266ea700 Mon Sep 17 00:00:00 2001 From: Chris Maffeo <cmaffeo2@illinois.edu> Date: Fri, 27 Sep 2019 13:17:47 -0500 Subject: [PATCH] Added sequence to vHelix parser --- mrdna/readers/__init__.py | 1 - mrdna/readers/polygon_mesh.py | 21 +++++++++++++++++++++ 2 files changed, 21 insertions(+), 1 deletion(-) diff --git a/mrdna/readers/__init__.py b/mrdna/readers/__init__.py index 31a6ced..4ac260c 100644 --- a/mrdna/readers/__init__.py +++ b/mrdna/readers/__init__.py @@ -13,7 +13,6 @@ def read_vhelix(maya_file, **model_parameters): maya_bases = parse_maya_file(maya_file) model = convert_maya_bases_to_segment_model( maya_bases, **model_parameters ) - model.set_sequence(m13seq*10) return model def read_list(infile,**model_parameters): diff --git a/mrdna/readers/polygon_mesh.py b/mrdna/readers/polygon_mesh.py index 2eb0015..7f96df1 100644 --- a/mrdna/readers/polygon_mesh.py +++ b/mrdna/readers/polygon_mesh.py @@ -13,6 +13,8 @@ maya_obj_registry_short = dict() class MayaObj(): """ Class representing a node in a Maya ascii file """ + seq_mapping = {str(i):c for c,i in zip('ATCG',range(4))} + def __init__(self, maya_lines): self.parent = None self.children = dict() @@ -59,6 +61,11 @@ class MayaObj(): Ry = rotationAboutAxis( [0,1,0], ay, normalizeAxis=False ) Rz = rotationAboutAxis( [0,0,1], az, normalizeAxis=False ) self.orientation = Rz.dot(Ry.dot(Rx)) + elif vals[1] == '".lb"': + try: + self.sequence = MayaObj.seq_mapping[vals[2][0]] + except: + pass def _strip_maya_string(string): return re.match('"\|?([a-zA-Z_0-9|]+)";?$', string).group(1) @@ -138,6 +145,7 @@ def ParseMayaConnection(line, base_dict): class MayaBase(MayaObj): def __init__(self, maya_lines): + self.sequence = '?' MayaObj.__init__(self, maya_lines) self._parse_maya_lines(maya_lines) self.basepair = None @@ -371,6 +379,10 @@ def convert_maya_bases_to_segment_model(maya_bases, **model_parameters): coord = np.zeros((len(maya_bases),3)) orientation = np.zeros((len(maya_bases),3,3)) bp,end3 = [-np.ones((len(maya_bases)),dtype=np.int) for i in range(2)] + sequence = [] + + seq_map = {'A':'T', 'C':'G', 'T':'A', 'G':'C'} + for i,b in zip(range(len(maya_bases)),maya_bases): coord[i] = b.get_position() @@ -380,6 +392,15 @@ def convert_maya_bases_to_segment_model(maya_bases, **model_parameters): if b.end3 is not None: end3[i] = base_to_index[ b.end3 ] + try: + s = b.sequence + if s == '?': + s = seq_map[s.basepair.sequence] + except: + s = 'T' + sequence.append( s ) + if 'sequence' not in model_parameters: + model_parameters['sequence'] = sequence return model_from_basepair_stack_3prime(coord, bp, None, end3, orientation=orientation, **model_parameters) -- GitLab