-
Notifications
You must be signed in to change notification settings - Fork 1
Expand file tree
/
Copy pathlcs.py
More file actions
60 lines (54 loc) · 1.91 KB
/
Copy pathlcs.py
File metadata and controls
60 lines (54 loc) · 1.91 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
#!/usr/bin/env python
# find an LCS (Longest Common Subsequence).
# *public domain*
def find_lcs_len(s1, s2):
m = [ [ 0 for x in s2 ] for y in s1 ]
for p1 in range(len(s1)):
for p2 in range(len(s2)):
if s1[p1] == s2[p2]:
if p1 == 0 or p2 == 0:
m[p1][p2] = 1
else:
m[p1][p2] = m[p1-1][p2-1]+1
elif m[p1-1][p2] < m[p1][p2-1]:
m[p1][p2] = m[p1][p2-1]
else: # m[p1][p2-1] < m[p1-1][p2]
m[p1][p2] = m[p1-1][p2]
return m[-1][-1]
def find_lcs(s1, s2):
# length table: every element is set to zero.
m = [ [ 0 for x in s2 ] for y in s1 ]
# direction table: 1st bit for p1, 2nd bit for p2.
d = [ [ None for x in s2 ] for y in s1 ]
# we don't have to care about the boundery check.
# a negative index always gives an intact zero.
for p1 in range(len(s1)):
for p2 in range(len(s2)):
if s1[p1] == s2[p2]:
if p1 == 0 or p2 == 0:
m[p1][p2] = 1
else:
m[p1][p2] = m[p1-1][p2-1]+1
d[p1][p2] = 3 # 11: decr. p1 and p2
elif m[p1-1][p2] < m[p1][p2-1]:
m[p1][p2] = m[p1][p2-1]
d[p1][p2] = 2 # 10: decr. p2 only
else: # m[p1][p2-1] < m[p1-1][p2]
m[p1][p2] = m[p1-1][p2]
d[p1][p2] = 1 # 01: decr. p1 only
(p1, p2) = (len(s1)-1, len(s2)-1)
# now we traverse the table in reverse order.
s = []
while (p1 >= 0 and p2 >= 0 and m[p1][p2]):
c = d[p1][p2]
if c == 3: s.append(s1[p1])
if c & 2: p2 -= 1
if c & 1: p1 -= 1
s.reverse()
return ''.join(s)
if __name__ == '__main__':
print find_lcs('BDCABA','ABCBDAB')
print find_lcs('XMJYAUZ','MZJAWXU')
print find_lcs('AAACCGTGAGTTATTCGTTCTAGAA','CACCCCTAAGGTACCTTTGGTTC')
#print find_lcs_len('XMJYAUZ','MZJAWXU')
#print find_lcs_len('AAACCGTGAGTTATTCGTTCTAGAA','CACCCCTAAGGTACCTTTGGTTC')