diff --git a/pysam/libcalignedsegment.pyx b/pysam/libcalignedsegment.pyx index 3071f375..20932293 100644 --- a/pysam/libcalignedsegment.pyx +++ b/pysam/libcalignedsegment.pyx @@ -75,6 +75,7 @@ from pysam.libcutils cimport force_bytes, force_str, \ charptr_to_str, charptr_to_bytes from pysam.libcutils cimport qualities_to_qualitystring, qualitystring_to_array, \ array_to_qualitystring +from libc.stdio cimport snprintf # Constants for binary tag conversion cdef char * htslib_types = 'cCsSiIf' @@ -1279,13 +1280,41 @@ cdef class AlignedSegment: empty string. ''' def __get__(self): - c = self.cigartuples - if c is None: + cdef uint32_t *cigar_p + cdef bam1_t *src + cdef uint32_t op, l + cdef int k + cdef int ret + cdef int pos = 0 + cdef int size = 16 + cdef char buf[16] + cdef char *s = buf + + src = self._delegate + if pysam_get_n_cigar(src) == 0: return None - # reverse order - else: - return "".join([ "%i%c" % (y,CODE2CIGAR[x]) for x,y in c]) + cigar_p = pysam_bam_get_cigar(src) + while True: + for k from 0 <= k < pysam_get_n_cigar(src): + op = cigar_p[k] & BAM_CIGAR_MASK + l = cigar_p[k] >> BAM_CIGAR_SHIFT + ret = snprintf(s + pos, size - pos, "%u%c", l, CODE2CIGAR[op]) + if ret >= size - pos: + pos = 0 + size = size * 2 + if s != buf: + free(s) + s = malloc(size) + break + pos += ret + else: + break + try: + return s[:pos].decode("utf8") + finally: + if s != buf: + free(s) def __set__(self, cigar): if cigar is None or len(cigar) == 0: self.cigartuples = [] diff --git a/tests/AlignedSegment_test.py b/tests/AlignedSegment_test.py index ea7886fa..8860ad77 100644 --- a/tests/AlignedSegment_test.py +++ b/tests/AlignedSegment_test.py @@ -1028,6 +1028,16 @@ def testStats(self): expected[1][i] = 1 self.assertEqual([list(x) for x in a.get_cigar_stats()], expected) + for i in range(1, 100): + cigarstring = "".join("10{}".format(x) + for x in iter("MIDNSHP=X")) * i + a.cigarstring = cigarstring + self.assertEqual(a.cigarstring, cigarstring) + expected = [[i * 10 for j in range(len("MIDNSHP=X"))] + [0, 0], + [i for j in range(len("MIDNSHP=X"))] + [0, 0]] + obtained = [list(x) for x in a.get_cigar_stats()] + self.assertEqual(obtained, expected) + a.cigarstring = "10M" a.set_tag("NM", 5) self.assertEqual(