laydi/scripts/geneontology/go-distance/godist.c

#include <string.h>
#include <errno.h>
#include <stdio.h>
#include <stdlib.h>
#include <search.h>
#include "godist.h"

void print_terms();
void add_link(char*, char*);
struct node* get_bp();

/* initialisation */
int godist_init() {
    /* Initialize hash table and array */
    hcreate(MAX_NODES);
    term_array_size = 0;
    link_count = 0;

    /* Read ontology terms from file */
    printf("Reading GO terms from go-terms.txt...");
    FILE *term_fd = fopen("go-terms.txt", "r");

    if (term_fd == NULL) {
        printf("cannot open file: go-terms.txt\n");
        exit(errno);
    }

    int i;
    while((i = godist_read_term(term_fd)) == 13) {
 /*       printf("%d\n", i);*/
    }
    fclose(term_fd);
    printf(" %d terms\n", term_array_size);

    /* Read ontology structure from file */
    printf("Reading GO structure from go-tree.txt...");
    FILE *tree_fd = fopen("go-tree.txt", "r");

    if (tree_fd == NULL) {
        printf("cannot open file: go-tree.txt\n");
        exit(errno);
    }

    while((i = godist_read_assoc(tree_fd)) == 2) {
	link_count++;
    }
    fclose(tree_fd);
    printf(" %d edges\n", link_count);

    accumulate_evidence(get_bp());
}

void godist_exit() {
    int i;
    for (i=0; i<term_array_size; i++) {
	free(term_array[i]);
    }
}

int godist_read_assoc(FILE *fd) {
    char term1[11], term2[11];
    int retval;
    retval = fscanf(fd, " %10s %10s ", term1, term2);
    if (retval != EOF) {
	add_link(term1, term2);
    }
    return retval;
}

int godist_read_term(FILE *fd) {
    char term[11];
    int ev[12];
    int i;
    ENTRY e, *res;
    int nread = fscanf(fd, " %10s %d %d %d %d %d %d %d %d %d %d %d %d ",
                       term, &ev[0], &ev[1], &ev[2], &ev[3], &ev[4], &ev[5],
                       &ev[6], &ev[7], &ev[8], &ev[9], &ev[10], &ev[11]);
    if (errno != 0) {
        printf("errno: %d\n", errno);
    }
    if (nread == 13) {
        struct node *n = (struct node*) malloc(sizeof(struct node));
	n->parentc = 0;
	n->childrenc = 0;
	for (i=0; i<12; i++) {
	    n->evidence[i] = ev[i];
	    n->acc_evidence[i] = 0;
	}
	strcpy(n->term, term);

	/* add to hash table */
	e.key = n->term;
	res = hsearch(e, ENTER);

        term_array[term_array_size++] = n;
    } 

    return nread;
}

/* distance functions */
float go_distance(char *term1, char *term2) {
    return 0.0;
}

void add_node(char *term, int parentc, struct node **parents) {
    
}

void add_link(char *parent_id, char *child_id) {
    ENTRY *ep, e;
    struct node *parent, *child;

    char key[11];
    strcpy(key, parent_id);
    e.key = key;
    ep = hsearch(e, FIND);
    if (!ep) {
	printf("Cannot find term %s\n", e.key);
	return;
    }
    parent = (struct node*) ep->key;

    strcpy(key, child_id);
    e.key = key;
    ep = hsearch(e, FIND);
    if (!ep) {
	printf("Cannot find term %s\n", e.key);
	return;
    }
    child = (struct node*) ep->key;

    if (parent->childrenc +1 > MAX_CHILDREN) {
        printf("FIXME: increase child count");
	return;
    }
    parent->children[parent->childrenc++] = child;
    child->parents[child->parentc++] = parent;
}

struct node *get_bp() {
    ENTRY e, *ep;
    e.key = "GO:0008150";
    ep = hsearch(e, FIND);
    if (ep)
	return ep;
    if (ep == ESRCH)
	printf("ESRCH\n");
    if (ep == ESRCH)
	printf("ENOMEM\n");
    return NULL;
}

void accumulate_evidence(struct node *n) {
    int i, j;
    for (i=0; i<(n->childrenc); i++) {
	accumulate_evidence(n->children[i]);
	for (j=0; j<12; j++)
            n->evidence[j] += n->children[i]->evidence[j];
    }
}

void print_terms() {
    int i;
    for (i=0; i<term_array_size; i++) {
	printf("%s\n", term_array[i]->term);
    }
}