Initial commit

MPI_packing.c

@@ -0,0 +1,173 @@
+#include "bison.h"
+
+/******************************************************************************
+*
+*   Take a BAM header and pack it into a single contiguous memory block. Store
+*   the resulting block and its size in an MPI_Header structure.
+*
+*   THE RESULT MUST BE free()d
+*
+*   bam_header_t *header: The header to store
+*
+*******************************************************************************/
+MPI_Header * pack_header(bam_header_t *header) {
+    size_t size = sizeof(int32_t); //n_targets
+    int32_t *pint32_t;
+    uint32_t *puint32_t;
+    char *pchar;
+    int *pint;
+    int i;
+    void *p;
+    MPI_Header *output = malloc(sizeof(MPI_Header));
+
+    //target_name
+    for(i=0; i<header->n_targets; i++) {
+        size += (sizeof(char) * (1+strlen(header->target_name[i])));
+    }
+
+    //target_len
+    size += sizeof(uint32_t) * header->n_targets;
+
+    //l_text
+    size += sizeof(int);
+
+    //text
+    size += sizeof(char) * (1 + header->l_text);
+
+    //Start copying, layout is n_targets,target_name[s],target_len[s],l_text,text
+    output->size = (int) size;
+    output->packed = malloc(size);
+    p = output->packed;
+
+    //n_targets
+    memcpy(p, (void *) &(header->n_targets), sizeof(int32_t));
+    pint32_t = (int32_t *) p;
+    p = (void *) (++pint32_t);
+
+    //target_name
+    for(i=0; i<header->n_targets; i++) {
+        memcpy(p, (void *) header->target_name[i], sizeof(char) * (1 + strlen(header->target_name[i])));
+        pchar = (char *) p;
+        p = (void *) (pchar+1+strlen(header->target_name[i]));
+    }
+    //target_len
+    memcpy(p, (void *) header->target_len, sizeof(uint32_t)*(header->n_targets));
+    puint32_t = (uint32_t *) p;
+    p = (void *) (puint32_t + header->n_targets);
+
+    //l_text
+    memcpy(p, (void *) &(header->l_text), sizeof(int));
+    pint = (int *) p;
+    p = (void *) ++pint;
+
+    //text
+    memcpy(p, (void *) (header->text), sizeof(char) * (1 + header->l_text));
+
+    return output;
+}
+
+/******************************************************************************
+*
+*   Unpack a header packed into an initialized bam_header_t
+*
+*   bam_header_t *header: The header to unpack into
+*   void *packed: The packed header
+*
+*******************************************************************************/
+void unpack_header(bam_header_t *header, void *packed) {
+    void *p = packed;
+    int i;
+    int *pint;
+    int32_t *pint32_t;
+    uint32_t *puint32_t;
+    char *pchar;
+    size_t strlength;
+
+    //n_targets
+    header->n_targets = *((int32_t *) packed);
+    pint32_t = (int32_t *) p;
+    p = (void *) (++pint32_t);
+
+    //**target_name
+    header->target_name = (char **) malloc(sizeof(char *) * (header->n_targets));
+    for(i=0; i<header->n_targets; i++) {
+        strlength = strlen((char *) p)+1;
+        header->target_name[i] = malloc(sizeof(char) * strlength);
+        memcpy((void *) (header->target_name[i]), p, sizeof(char)*strlength);
+        pchar = (char *) p;
+        p = (void *) (pchar+strlength);
+    }
+
+    //target_len
+    header->target_len = malloc(sizeof(uint32_t) * (header->n_targets));
+    for(i=0; i<header->n_targets; i++) {
+        header->target_len[i] = *((uint32_t *) p);
+        puint32_t = (uint32_t *) p;
+        p = (void *) ++puint32_t;
+    }
+
+    //l_text
+    header->l_text = *((int *) p);
+    pint = (int *) p;
+    p = (void *) ++pint;
+
+    //text
+    header->text = (char *) malloc(sizeof(char) * (header->l_text+1));
+    memcpy((void *) (header->text), p, sizeof(char) * (header->l_text + 1));
+}
+
+/******************************************************************************
+*
+*   Take a BAM read and pack it into a single contiguous memory block. Store
+*   the resulting block and its size in an MPI_Read structure.
+*
+*   THE RESULT MUST BE free()d
+*
+*   bam1_t *read: The read to store
+*
+*******************************************************************************/
+MPI_read * pack_read(bam1_t *read, MPI_read *output) {
+    bam1_t *pbam1_t;
+    int needed_size, m_data = read->m_data;
+
+    needed_size = (int) (sizeof(bam1_t) + m_data);
+    if(output->size == 0) {
+        output->packed = malloc((size_t) needed_size);
+        output->size = needed_size;
+    } else if(needed_size > output->size) {
+        output->packed = realloc(output->packed, (size_t) needed_size);
+        output->size = needed_size;
+    }
+    memcpy((void *) output->packed, (void *) read, sizeof(bam1_t));
+    pbam1_t = output->packed;
+    pbam1_t++;
+    memcpy((void *) pbam1_t, (void *) read->data, m_data);
+    return output;
+}
+
+/******************************************************************************
+*
+*   Unpack a packed read into an initialized bam1_t read.
+*
+*   bam1_t *read: The read to unpack into
+*   void *packed: The packed read
+*
+*******************************************************************************/
+bam1_t *unpack_read(bam1_t *read, void *packed) {
+    bam1_t *pbam1_t = packed;
+    uint8_t *pdata = (uint8_t *) (pbam1_t+1);
+    uint8_t *newdata;
+
+    pbam1_t->data = pdata;
+    if(read != NULL) bam_destroy1(read);
+    read = bam_init1();
+    read->core = pbam1_t->core;
+    read->l_aux = pbam1_t->l_aux;
+    read->m_data = pbam1_t->m_data;
+    read->data_len= pbam1_t->data_len;
+    newdata = (uint8_t *) malloc(read->m_data);
+    memcpy((void *) newdata, (void *) pdata, read->m_data);
+    read->data = newdata;
+
+    return read;
+}

Makefile

@@ -0,0 +1,82 @@
+WORK=/home/ryand#This should be changed to match your needs
+PREFIX = $(WORK)/bin
+CC = mpicc
+INCLUDE_DIRS = -I$(WORK)/include #This should be were samtools was compiled -I/path/to/samtools/compilation
+LIB_DIRS = -L$(WORK)/lib #As above, but -L/path/to/samtools/compilation
+OPTS = -Wall -O3 #-DDEBUG #-DNOTHROTTLE -g
+MPI = -lmpich -lmpl #This is usually appropriate for mpich2
+#MPI = #This is appropriate for mvapich2
+#MPI = -lmpi #This is usually appropriate for openmpi
+
+#Don't edit below here unless you know what you're doing!
+
+OBJS = aux.o fastq.o genome.o slurp.o master.o common.o MPI_packing.o worker.o
+HERD_OBJS = herd/fastq.o herd/master.o herd/MPI_packing.o herd/slurp.o herd/worker.o herd/writer.o
+
+.SUFFIXES:.c .o
+
+all: align index extractor mbias markduplicates
+
+.c.o:
+	$(CC) -c $(OPTS) $(INCLUDE_DIRS) $< -o $@
+
+markduplicates:
+	$(CC) $(OPTS) $(INCLUDE_DIRS) $(LIB_DIRS) -o bison_markduplicates markduplicates.c -lpthread -lbam -lz
+
+mbias:
+	$(CC) $(OPTS) $(INCLUDE_DIRS) $(LIB_DIRS) -o bison_mbias mbias.c -lpthread -lbam -lz
+
+index:
+	$(CC) $(OPTS) -o bison_index index.c -lpthread
+
+align: $(OBJS)
+	$(CC) -c $(OPTS) $(INCLUDE_DIRS) main.c -o main.o
+	$(CC) $(OPTS) $(OBJS) main.o -o bison $(LIB_DIRS) -lm -lpthread $(MPI) -lbam -lz
+
+extractor:
+	$(CC) -c $(OPTS) $(INCLUDE_DIRS) common.c -o common.o
+	$(CC) -c $(OPTS) $(INCLUDE_DIRS) methylation_extractor.c -o methylation_extractor.o
+	$(CC) $(OPTS) $(LIB_DIRS) common.o methylation_extractor.o -o bison_methylation_extractor -lpthread -lbam -lz
+
+#Don't compile herd by default
+herd:  $(OBJS) $(HERD_OBJS)
+	$(CC) -c $(OPTS) $(INCLUDE_DIRS) herd/main.c -o herd/main.o
+	$(CC) $(OPTS) $(OBJS) $(HERD_OBJS) herd/main.o -o bison_herd $(LIB_DIRS) -lm -lpthread $(MPI) -lbam -lz
+
+#Auxiliary programs, don't compile by default
+auxiliary:	merge_CpGs bedGraph2methylKit make_reduced_genome aux_python_scripts CpG_coverage
+
+aux_python_scripts:
+	cp -f auxiliary/bedGraph2BSseq.py ./
+	cp -f auxiliary/merge_bedGraphs.py ./
+
+CpG_coverage:	common.o
+	$(CC) -c $(OPTS) $(INCLUDE_DIRS) auxiliary/CpG_coverage.c -o auxiliary/CpG_coverage.o
+	$(CC) $(OPTS) $(LIB_DIRS) common.o auxiliary/CpG_coverage.o -o bison_CpG_coverage
+
+merge_CpGs:	common.o
+	$(CC) -c $(OPTS) $(INCLUDE_DIRS) auxiliary/merge_CpGs.c -o auxiliary/merge_CpGs.o
+	$(CC) $(OPTS) $(LIB_DIRS) common.o auxiliary/merge_CpGs.o -o bison_merge_CpGs
+
+bedGraph2methylKit:common.o
+	$(CC) -c $(OPTS) $(INCLUDE_DIRS) auxiliary/bedGraph2methylKit.c -o auxiliary/bedGraph2methylKit.o
+	$(CC) $(OPTS) $(LIB_DIRS) common.o auxiliary/bedGraph2methylKit.o -o bedGraph2methylKit
+
+make_reduced_genome:
+	$(CC) $(OPTS) $(LIB_DIRS) auxiliary/make_reduced_genome.c -o make_reduced_genome
+
+install :
+	mv bison_* $(PREFIX)/ 
+	chmod a+x Rscripts/*
+	cp Rscripts/* $(PREFIX)/
+	if [ -f bison ]; then mv bison $(PREFIX)/ ; fi;
+	if [ -f bedGraph2methylKit ]; then mv bedGraph2methylKit $(PREFIX)/ ; fi;
+	if [ -f bedGraph2BSseq.py ]; then chmod a+x bedGraph2BSseq.py ; mv bedGraph2BSseq.py $(PREFIX)/ ; fi;
+	if [ -f merge_bedGraphs.py ]; then chmod a+x merge_bedGraphs.py ; mv merge_bedGraphs.py $(PREFIX)/ ; fi;
+	if [ -f check_accuracy ]; then mv check_accuracy $(PREFIX)/ ; fi;
+	if [ -f make_reduced_genome ]; then mv make_reduced_genome $(PREFIX)/ ; fi;
+
+clean:
+	rm -f *.o bison bison_* bedGraph2methylKit check_accuracy make_reduced_genome bedGraph2BSseq.py
+	rm -f herd/*.o
+	rm -f auxiliary/*.o