lcs.c File Reference

Longest Common Subsequence algorithm More...

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <assert.h>

Include dependency graph for lcs.c:

Enumerations
enum	{ LEFT , UP , DIAG }

Functions
void	lcslen (const char *s1, const char *s2, int l1, int l2, int **L, int **B)
	Computes LCS between s1 and s2 using a dynamic-programming approach.
char *	lcsbuild (const char *s1, int l1, int l2, int **L, int **B)
	Builds the LCS according to B using a traceback approach.
static void	test ()
	Self-test implementations.
int	main (int argc, char *argv[])
	Main function.

Detailed Description

Longest Common Subsequence algorithm

From Wikipedia: The longest common subsequence (LCS) problem is the problem of finding the longest subsequence common to all sequences in a set of sequences (often just two sequences).

Author

Kurtz

Enumeration Type Documentation

anonymous enum

anonymous enum

{LEFT, UP, DIAG};

Function Documentation

lcsbuild()

char * lcsbuild	(	const char *	s1,
		int	l1,
		int	l2,
		int **	L,
		int **	B
	)

Builds the LCS according to B using a traceback approach.

Parameters

s1	first null-terminated string
l1	length of s1
l2	length of s2
L	matrix of size l1 x l2
B	matrix of size l1 x l2

Returns

lcs longest common subsequence

                                                                 {
    int  i, j, lcsl;
    char    *lcs;
    lcsl = L[l1][l2];
    
    /* my lcs is at least the empty symbol */
    lcs = (char *)calloc(lcsl+1, sizeof(char)); /* null-terminated \0 */
    if (!lcs) {
        perror("calloc: ");
        return NULL;
    }
 
    i = l1, j = l2;
    while (i > 0 && j > 0) {
        /* walk the matrix backwards */
        if (B[i][j] == DIAG) {
            lcs[--lcsl] = s1[i-1];
            i = i - 1;
            j = j - 1;
        }
        else if (B[i][j] == LEFT)
        {
            j = j - 1;
        }
        else
        {
            i = i - 1;
        }
    }
    return lcs;
}

lcslen()

void lcslen	(	const char *	s1,
		const char *	s2,
		int	l1,
		int	l2,
		int **	L,
		int **	B
	)

Computes LCS between s1 and s2 using a dynamic-programming approach.

Parameters

s1	first null-terminated string
s2	second null-terminated string
l1	length of s1
l2	length of s2
L	matrix of size l1 x l2
B	matrix of size l1 x l2

Returns

void

                                                                              {
    /* B is the directions matrix
       L is the LCS matrix */
    int i, j;
 
    /* loop over the simbols in my sequences
       save the directions according to the LCS */
    for (i = 1; i <= l1; ++i) {
        for (j = 1; j <= l2; ++j) {
            if (s1[i-1] == s2[j-1]) {
                L[i][j] = 1 + L[i-1][j-1];
                B[i][j] = DIAG;
            }
            else if (L[i-1][j] < L[i][j-1]) {
                L[i][j] = L[i][j-1];
                B[i][j] = LEFT;
            }
            else {
                L[i][j] = L[i-1][j];
                B[i][j] = UP;
            }
        }
    }
}

main()

int main	(	int	argc,
		char *	argv[]
	)

Main function.

Parameters

argc	commandline argument count (ignored)
argv	commandline array of arguments (ignored)

Returns

0 on exit

                                 {
    test();  // run self-test implementations
    return 0;
}

Here is the call graph for this function:

test()

static void test ( void  )

static

Self-test implementations.

Returns

void

                   {
    /* https://en.wikipedia.org/wiki/Subsequence#Applications */
    int **L, **B, j, l1, l2;
    
    char *s1 = "ACGGTGTCGTGCTATGCTGATGCTGACTTATATGCTA";
    char *s2 = "CGTTCGGCTATCGTACGTTCTATTCTATGATTTCTAA";
    char *lcs;
    
    l1 = strlen(s1);
    l2 = strlen(s2);
    
    L = (int **)calloc(l1+1, sizeof(int *));
    B = (int **)calloc(l1+1, sizeof(int *));
    
    if (!L) {
        perror("calloc: ");
        exit(1);
    }
    if (!B) {
        perror("calloc: ");
        exit(1);
    }
    for (j = 0; j <= l1; j++) {
        L[j] = (int *)calloc(l2+1, sizeof(int));
        if (!L[j]) {
            perror("calloc: ");
            exit(1);
        }
        B[j] = (int *)calloc(l2+1, sizeof(int));
        if (!L[j]) {
            perror("calloc: ");
            exit(1);
        }
    }
    
    lcslen(s1, s2, l1, l2, L, B);
    lcs = lcsbuild(s1, l1, l2, L, B);
    
    assert(L[l1][l2] == 27);
    assert(strcmp(lcs, "CGTTCGGCTATGCTTCTACTTATTCTA") == 0);
    
    printf("S1: %s\tS2: %s\n", s1, s2);
    printf("LCS len:%3d\n", L[l1][l2]);
    printf("LCS: %s\n", lcs);
 
    free(lcs);
    for (j = 0; j <= l1; j++)
    {
        free(L[j]), free(B[j]);
    }
    free(L);
    free(B);
 
    printf("All tests have successfully passed!\n");
}

Here is the call graph for this function: