Commits

Anonymous committed 9107c4c

Corrected problem in maf_extract_ranges_indexed, which failed when a maf block contains more than one line with the desired reference sequence. The failure manifested as an abort, from the range error test in coord_to_col() in align/core.py. The correction is two-fold. First, we add a test in the script co that coord_to_col is not called with an empty interval. Second, we change the script to loop over *all* interesting reference intervals, rather than the first.

Comments (0)

Files changed (2)

             if c.src == src: return c
         return None
 
+    def get_components_by_src( self, src ):
+        for c in self.components:
+            if c.src == src: yield c
+
     def get_component_by_src_start( self, src ):
         for c in self.components:
             if c.src.startswith( src ): return c

scripts/maf_extract_ranges_indexed.py

         # Write each intersecting block
         if chop:
             for block in blocks: 
-                ref = block.get_component_by_src( src )
-                slice_start = max( start, ref.get_forward_strand_start() )
-                slice_end = min( end, ref.get_forward_strand_end() )
-                sliced = block.slice_by_component( ref, slice_start, slice_end ) 
-                # If the block is shorter than the minimum allowed size, stop
-                if mincols and ( sliced.text_size < mincols ):
-                    continue
-                # If the reference component is empty, don't write the block
-                if sliced.get_component_by_src( src ).size < 1:
-                    continue
-                # Keep only components that are not empty
-                sliced.components = [ c for c in sliced.components if c.size > 0 ]
-                # Reverse complement if needed
-                if ( strand != None ) and ( ref.strand != strand ): 
-                    sliced = sliced.reverse_complement()
-                # Write the block
-                out.write( sliced )
+                for ref in block.get_components_by_src( src ):
+                    slice_start = max( start, ref.get_forward_strand_start() )
+                    slice_end = min( end, ref.get_forward_strand_end() )
+                    if (slice_end <= slice_start): continue
+                    sliced = block.slice_by_component( ref, slice_start, slice_end ) 
+                    # If the block is shorter than the minimum allowed size, stop
+                    if mincols and ( sliced.text_size < mincols ):
+                        continue
+                    # If the reference component is empty, don't write the block
+                    if sliced.get_component_by_src( src ).size < 1:
+                        continue
+                    # Keep only components that are not empty
+                    sliced.components = [ c for c in sliced.components if c.size > 0 ]
+                    # Reverse complement if needed
+                    if ( strand != None ) and ( ref.strand != strand ): 
+                        sliced = sliced.reverse_complement()
+                    # Write the block
+                    out.write( sliced )
         else:
             for block in blocks:
                 out.write( block )