Commit 0758865d by Eric Coissac

Fix #45

parent 3ccff024
......@@ -117,10 +117,10 @@ cpdef dict codon(char aa, bint direct=True):
if direct:
# value = <object>PyDict_GetItemString(UniversalGeneticCode,letter)
value = UniversalGeneticCode[letter]
value = UniversalGeneticCode.get(letter,())
else:
# value = <object>PyDict_GetItemString(CompUniversalGeneticCode,letter)
value = CompUniversalGeneticCode[letter]
value = CompUniversalGeneticCode.get(letter,())
return listCodons2Dict(value)
......@@ -838,14 +838,15 @@ cdef class AhoCorasick:
cdef class ProtAhoCorasick(AhoCorasick):
cpdef int addSequence(self,bytes sequence, object seqid=1, size_t kup=4):
cdef int position
cdef int k
cdef char* csequence
cdef char letter[2]
cdef dict c
cdef int osize
cdef int position
cdef int k
cdef char* csequence
cdef char letter[2]
cdef dict c
cdef dict cn
cdef int osize
cdef size_t erreur
cdef int id
cdef int id
if self._finalized:
return 0
......@@ -874,15 +875,23 @@ cdef class ProtAhoCorasick(AhoCorasick):
c = codon(csequence[position+kup-1])
for k in range(position+kup-2,position-1,-1):
c = bindCodons(codon(csequence[k]),c)
c = word2automata([w for w in enumerateword(c) if not homopolymer(w)])
erreur = self.addAutomata(c, id, (position+kup) * 3 - 1)
if erreur==0:
return 0
if c:
erreur=0
for k in range(position+kup-2,position-1,-1):
cn = codon(csequence[k])
if cn:
c = bindCodons(cn,c)
else:
erreur=1
if erreur==0:
c = word2automata([w for w in enumerateword(c) if not homopolymer(w)])
erreur = self.addAutomata(c, id, (position+kup) * 3 - 1)
if erreur==0:
return 0
# Backtranslate according to the reverse strand
for position in range(len(sequence)-1,kup-2,-1):
......@@ -891,14 +900,22 @@ cdef class ProtAhoCorasick(AhoCorasick):
# anticodon variants
c = codon(csequence[position-kup+1],False)
for k in range(position-kup+2,position+1):
c = bindCodons(codon(csequence[k],False),c)
c = word2automata([w for w in enumerateword(c) if not homopolymer(w)])
erreur = self.addAutomata(c, -id, position * 3 + 2)
if erreur==0:
return 0
if c:
erreur=0
for k in range(position-kup+2,position+1):
cn = codon(csequence[k],False)
if cn:
c = bindCodons(cn,c)
else:
erreur=1
if erreur==0:
c = word2automata([w for w in enumerateword(c) if not homopolymer(w)])
erreur = self.addAutomata(c, -id, position * 3 + 2)
if erreur==0:
return 0
return self._statecount - osize
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment