Commit 5c917133 by Eric Coissac

Patch a bug in ecoTag

parent da92eca1
...@@ -4,7 +4,7 @@ from obitools.fasta import fastaNucIterator,formatFasta ...@@ -4,7 +4,7 @@ from obitools.fasta import fastaNucIterator,formatFasta
#from obitools.align.ssearch import ssearchIterator #from obitools.align.ssearch import ssearchIterator
from obitools.utils.bioseq import uniqSequence,sortSequence from obitools.utils.bioseq import uniqSequence,sortSequence
from obitools.align import lenlcs from obitools.align import lenlcs,ALILEN
from obitools.options.taxonomyfilter import addTaxonomyDBOptions,loadTaxonomyDatabase from obitools.options.taxonomyfilter import addTaxonomyDBOptions,loadTaxonomyDatabase
from obitools.ecobarcode.databases import referenceDBIterator from obitools.ecobarcode.databases import referenceDBIterator
...@@ -123,13 +123,8 @@ def lcsIterator(entries,db,options): ...@@ -123,13 +123,8 @@ def lcsIterator(entries,db,options):
maxid = (None,0.0) maxid = (None,0.0)
minid = options.minimum minid = options.minimum
for d in db: for d in db:
lcs,lali = lenlcs(seq,d,minid,normalized=True) lcs,lali = lenlcs(seq,d,minid,normalized=True,reference=ALILEN)
# lcs,lali = lenlcs(seq,d,minid,normalized=True, large=options.large)
if lcs > maxid[1]: if lcs > maxid[1]:
# try:
# print "@@@@@@@@@@@>>>>>>",(maxid[0].id,maxid[1]),(d.id,lcs)
# except:
# print "@@@@@@@@@@@>>>>>>",(None,maxid[1]),(d.id,lcs)
maxid = (d,lcs) maxid = (d,lcs)
minid = maxid[1] ** options.shape minid = maxid[1] ** options.shape
results.append((d,lcs)) results.append((d,lcs))
...@@ -145,7 +140,7 @@ def lcsIteratorSelf(entries,db,options): ...@@ -145,7 +140,7 @@ def lcsIteratorSelf(entries,db,options):
maxid = (None,0.0) maxid = (None,0.0)
minid = options.minimum minid = options.minimum
for d in db: for d in db:
lcs,lali = lenlcs(seq,d,minid,normalized=True) lcs,lali = lenlcs(seq,d,minid,normalized=True,reference=ALILEN)
if lcs > maxid[1]: if lcs > maxid[1]:
maxid = (d,lcs) maxid = (d,lcs)
minid = maxid[1] minid = maxid[1]
...@@ -154,7 +149,7 @@ def lcsIteratorSelf(entries,db,options): ...@@ -154,7 +149,7 @@ def lcsIteratorSelf(entries,db,options):
results.append(maxid) results.append(maxid)
for d in db: for d in db:
if d.id != maxid[0].id: if d.id != maxid[0].id:
lcs,lali = lenlcs(maxid[0],d,maxid[1],normalized=True) lcs,lali = lenlcs(maxid[0],d,maxid[1],normalized=True,reference=ALILEN)
if lcs >= maxid[1]: if lcs >= maxid[1]:
results.append((d,lcs)) results.append((d,lcs))
......
...@@ -171,7 +171,7 @@ if __name__=='__main__': ...@@ -171,7 +171,7 @@ if __name__=='__main__':
taxonok+=1 taxonok+=1
if taxoncount: if taxoncount:
print '%-20s\t%10d\t%10d\t%8.2f' % (rank,taxonok,taxoncount,float(taxonok)/taxoncount*100) print '%-20s\t%10d\t%10d\t%8.2f' % (rank,taxonok,taxoncount,float(taxonok)/taxoncount*100)
# if len (lostSeqs) > 0: # if len (lostSeqs) > 0:
# print "Lost Sequences:" # print "Lost Sequences:"
# print lostSeqs # print lostSeqs
......
...@@ -2,7 +2,7 @@ ...@@ -2,7 +2,7 @@
from _nws import NWS from _nws import NWS
from _upperbond import indexSequences from _upperbond import indexSequences
from _lcs import LCS,lenlcs from _lcs import LCS,lenlcs,ALILEN,MAXLEN,MINLEN
from _assemble import DirectAssemble, ReverseAssemble from _assemble import DirectAssemble, ReverseAssemble
from _qsassemble import QSolexaDirectAssemble,QSolexaReverseAssemble from _qsassemble import QSolexaDirectAssemble,QSolexaReverseAssemble
from _rassemble import RightDirectAssemble as RightReverseAssemble from _rassemble import RightDirectAssemble as RightReverseAssemble
......
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