composition.c

#include <tcl.h>

#include "IO.h"                 /* GapIO */
#include "gap_globals.h"        /* consensus/quality cutoffs */
#include "qual.h"               /* calc_consensus() */
#include "cli_arg.h"            /* cli_arg, parse_args() */

static int tcl_composition(ClientData clientData, Tcl_Interp *interp,
                           int argc, char **argv);
static char *doit(GapIO *io, int contig, int lreg, int rreg);

/*
 * This is called when the library is dynamically linked in with the calling
 * program. Use it to initialise any tables and to register the necessary
 * commands.
 */
int Composition_Init(Tcl_Interp *interp) {
    if (NULL == Tcl_CreateCommand(interp,
                                  "composition",
                                  tcl_composition,
                                  (ClientData) NULL,
                                  (Tcl_CmdDeleteProc *) NULL))
        return TCL_ERROR;

    return TCL_OK;
}

/*
 * The composition itself.
 * This is called with an argc and argv in much the same way that main()
 * is. We can either parse them ourselves, our use the gap parse_args
 * utility routine.
 */
static int tcl_composition(ClientData clientData, Tcl_Interp *interp,
                           int argc, char **argv) {
    int num_contigs;
    contig_list_t *contigs = NULL;
    char *result;
    int i;
    Tcl_DString dstr;

    /* A structure definition to store the arguments in */
    typedef struct {
        GapIO *io;
        char *ident;
    } test_args;

    /* The mapping of the argument strings to our structure above */
    test_args args;
    cli_args a[] = {
        {"-io",       ARG_IO,  1, NULL, offsetof(test_args, io)},
        {"-contigs",  ARG_STR, 1, NULL, offsetof(test_args, ident)},
        {NULL,      0,       0, NULL, 0}
    };

    /*
     * First things first, add a header to the output window. This shows the
     * date and function name.
     */
    vfuncheader("test command");

    /* Parse the arguments */
    if (-1 == gap_parse_args(a, &args, argc, argv)) {
        return TCL_ERROR;
    }

    active_list_contigs(args.io, args.ident, &num_contigs, &contigs);
    if (num_contigs == 0) {
        xfree(contigs);
        return TCL_OK;
    }

    /* Do the actual work */
    Tcl_DStringInit(&dstr);
    for (i = 0; i < num_contigs; i++) {
        result = doit(args.io, contigs[i].contig, contigs[i].start,
                      contigs[i].end);
        if (NULL == result) {
            xfree(contigs);
            return TCL_ERROR;
        }

        Tcl_DStringAppendElement(&dstr, result);
    }

    Tcl_DStringResult(interp, &dstr);

    xfree(contigs);
    return TCL_OK;
}

/*
 * Our main work horse. For something to do as an example we'll output
 * the sequence composition of the contig in the given range.
 */
static char *doit(GapIO *io, int contig, int lreg, int rreg) {
    static char result[1024];
    char *consensus;
    int i, n[5];

    if (0 == lreg && 0 == rreg) {
        rreg = io_clength(io, contig);
        lreg = 1;
    }

    if (NULL == (consensus = (char *)xmalloc(rreg-lreg+1)))
        return NULL;

    if (-1 == calc_consensus(contig, lreg, rreg, CON_SUM,
                             consensus, NULL, NULL, NULL,
                             consensus_cutoff, quality_cutoff,
                             database_info, (void *)io)) {
        xfree(consensus);
        return NULL;
    }

    n[0] = n[1] = n[2] = n[3] = n[4] = 0;
    for (i = 0; i <= rreg - lreg; i++) {
        switch(consensus[i]) {
        case 'a':
        case 'A':
            n[0]++;
            break;

        case 'c':
        case 'C':
            n[1]++;
            break;

        case 'g':
        case 'G':
            n[2]++;
            break;

        case 't':
        case 'T':
            n[3]++;
            break;

        default:
            n[4]++;
        }
    }

    /* Return the information */
    sprintf(result, "%d %d %d %d %d %d",
            rreg - lreg + 1, n[0], n[1], n[2], n[3], n[4]);

    xfree(consensus);

    return result;
}