# Example for: sequence_db.search()

# This will search the MODELLER database of representative protein chains
# for chains similar to the specified sequence.

log.verbose()
env = environ()

# Read in the standard database included in the Modeller distribution
sdb = sequence_db(env, seq_database_file='$(LIB)/CHAINS_all.seq',
                  chains_list='very-short-for-test.cod')

# Read in the query sequence in alignment format
aln = alignment(env, file='toxin.ali', align_codes='1nbt')

sdb.search(aln, search_randomizations=20, # should use 100 in real life
           seq_database_file='$(LIB)/CHAINS_all.seq',
           off_diagonal=9999, gap_penalties_1d=(-800, -400),
           signif_cutoff=(1.5, 5.0))

aln.malign()
aln.write(file='toxin-search.pap', alignment_format='PAP')