--- /dev/null
+#include "pysam.h"
+
+/*\r
+\r
+bam_cat -- efficiently concatenates bam files\r
+\r
+bam_cat can be used to concatenate BAM files. Under special\r
+circumstances, it can be used as an alternative to 'samtools merge' to\r
+concatenate multiple sorted files into a single sorted file. For this\r
+to work each file must be sorted, and the sorted files must be given\r
+as command line arguments in order such that the final read in file i\r
+is less than or equal to the first read in file i+1.\r
+\r
+This code is derived from the bam_reheader function in samtools 0.1.8\r
+and modified to perform concatenation by Chris Saunders on behalf of\r
+Illumina.\r
+\r
+\r
+########## License:\r
+\r
+The MIT License\r
+\r
+Original SAMtools work copyright (c) 2008-2009 Genome Research Ltd.\r
+Modified SAMtools work copyright (c) 2010 Illumina, Inc.\r
+\r
+Permission is hereby granted, free of charge, to any person obtaining a copy\r
+of this software and associated documentation files (the "Software"), to deal\r
+in the Software without restriction, including without limitation the rights\r
+to use, copy, modify, merge, publish, distribute, sublicense, and/or sell\r
+copies of the Software, and to permit persons to whom the Software is\r
+furnished to do so, subject to the following conditions:\r
+\r
+The above copyright notice and this permission notice shall be included in\r
+all copies or substantial portions of the Software.\r
+\r
+THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR\r
+IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,\r
+FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE\r
+AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER\r
+LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,\r
+OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN\r
+THE SOFTWARE.\r
+\r
+*/\r
+\r
+\r
+/*\r
+makefile:\r
+"""\r
+CC=gcc\r
+CFLAGS+=-g -Wall -O2 -D_FILE_OFFSET_BITS=64 -D_USE_KNETFILE -I$(SAMTOOLS_DIR)\r
+LDFLAGS+=-L$(SAMTOOLS_DIR)\r
+LDLIBS+=-lbam -lz\r
+\r
+all:bam_cat\r
+"""\r
+*/\r
+\r
+\r
+#include <stdio.h>\r
+#include <stdlib.h>\r
+#include <unistd.h>\r
+\r
+#include "bgzf.h"\r
+#include "bam.h"\r
+\r
+#define BUF_SIZE 0x10000\r
+\r
+#define GZIPID1 31\r
+#define GZIPID2 139\r
+\r
+#define BGZF_EMPTY_BLOCK_SIZE 28\r
+\r
+\r
+int bam_cat(int nfn, char * const *fn, const bam_header_t *h, const char* outbam)\r
+{\r
+ BGZF *fp;\r
+ FILE* fp_file;\r
+ uint8_t *buf;\r
+ uint8_t ebuf[BGZF_EMPTY_BLOCK_SIZE];\r
+ const int es=BGZF_EMPTY_BLOCK_SIZE;\r
+ int i;\r
+ \r
+ fp = strcmp(outbam, "-")? bgzf_open(outbam, "w") : bgzf_fdopen(fileno(stdout), "w");\r
+ if (fp == 0) {\r
+ fprintf(pysamerr, "[%s] ERROR: fail to open output file '%s'.\n", __func__, outbam);\r
+ return 1;\r
+ }\r
+ if (h) bam_header_write(fp, h);\r
+ \r
+ buf = (uint8_t*) malloc(BUF_SIZE);\r
+ for(i = 0; i < nfn; ++i){\r
+ BGZF *in;\r
+ bam_header_t *old;\r
+ int len,j;\r
+ \r
+ in = strcmp(fn[i], "-")? bam_open(fn[i], "r") : bam_dopen(fileno(stdin), "r");\r
+ if (in == 0) {\r
+ fprintf(pysamerr, "[%s] ERROR: fail to open file '%s'.\n", __func__, fn[i]);\r
+ return -1;\r
+ }\r
+ if (in->open_mode != 'r') return -1;\r
+ \r
+ old = bam_header_read(in);\r
+ if (h == 0 && i == 0) bam_header_write(fp, old);\r
+ \r
+ if (in->block_offset < in->block_length) {\r
+ bgzf_write(fp, in->uncompressed_block + in->block_offset, in->block_length - in->block_offset);\r
+ bgzf_flush(fp);\r
+ }\r
+ \r
+ j=0;\r
+#ifdef _USE_KNETFILE\r
+ fp_file=fp->x.fpw;\r
+ while ((len = knet_read(in->x.fpr, buf, BUF_SIZE)) > 0) {\r
+#else \r
+ fp_file=fp->file;\r
+ while (!feof(in->file) && (len = fread(buf, 1, BUF_SIZE, in->file)) > 0) {\r
+#endif\r
+ if(len<es){\r
+ int diff=es-len;\r
+ if(j==0) {\r
+ fprintf(pysamerr, "[%s] ERROR: truncated file?: '%s'.\n", __func__, fn[i]);\r
+ return -1;\r
+ }\r
+ fwrite(ebuf, 1, len, fp_file);\r
+ memcpy(ebuf,ebuf+len,diff);\r
+ memcpy(ebuf+diff,buf,len);\r
+ } else {\r
+ if(j!=0) fwrite(ebuf, 1, es, fp_file);\r
+ len-= es;\r
+ memcpy(ebuf,buf+len,es);\r
+ fwrite(buf, 1, len, fp_file);\r
+ }\r
+ j=1;\r
+ }\r
+\r
+ /* check final gzip block */\r
+ {\r
+ const uint8_t gzip1=ebuf[0];\r
+ const uint8_t gzip2=ebuf[1];\r
+ const uint32_t isize=*((uint32_t*)(ebuf+es-4));\r
+ if(((gzip1!=GZIPID1) || (gzip2!=GZIPID2)) || (isize!=0)) {\r
+ fprintf(pysamerr, "[%s] WARNING: Unexpected block structure in file '%s'.", __func__, fn[i]);\r
+ fprintf(pysamerr, " Possible output corruption.\n");\r
+ fwrite(ebuf, 1, es, fp_file);\r
+ }\r
+ }\r
+ bam_header_destroy(old);\r
+ bgzf_close(in);\r
+ }\r
+ free(buf);\r
+ bgzf_close(fp);\r
+ return 0;\r
+}\r
+\r
+\r
+\r
+int main_cat(int argc, char *argv[])\r
+{\r
+ bam_header_t *h = 0;\r
+ char *outfn = 0;\r
+ int c, ret;\r
+ while ((c = getopt(argc, argv, "h:o:")) >= 0) {\r
+ switch (c) {\r
+ case 'h': {\r
+ tamFile fph = sam_open(optarg);\r
+ if (fph == 0) {\r
+ fprintf(pysamerr, "[%s] ERROR: fail to read the header from '%s'.\n", __func__, argv[1]);\r
+ return 1;\r
+ }\r
+ h = sam_header_read(fph);\r
+ sam_close(fph);\r
+ break;\r
+ }\r
+ case 'o': outfn = strdup(optarg); break;\r
+ }\r
+ }\r
+ if (argc - optind < 2) {\r
+ fprintf(pysamerr, "Usage: samtools cat [-h header.sam] [-o out.bam] <in1.bam> <in2.bam> [...]\n");\r
+ return 1;\r
+ }\r
+ ret = bam_cat(argc - optind, argv + optind, h, outfn? outfn : "-");\r
+ free(outfn);\r
+ return ret;\r
+}\r
projects / pysam.git / blobdiff