cb9df600fd62bdbb03fc0f1bee61a6f2b0440e09 aamp Sun Jul 17 09:56:30 2011 -0700 made a spinoff function of bamFetch() that also handles an already-open bam file/index for cases where fetch is done within a loop diff --git src/lib/bamFile.c src/lib/bamFile.c index 0f9a9e3..1ef3588 100644 --- src/lib/bamFile.c +++ src/lib/bamFile.c @@ -81,71 +81,81 @@ errAbort("Failed to open %s%s", fileOrUrl, urlWarning->string); } return fh; } void bamClose(samfile_t **pSamFile) /* Close down a samefile_t */ { if (pSamFile != NULL) { samclose(*pSamFile); *pSamFile = NULL; } } +void bamFetchAlreadyOpen(samfile_t *samfile, bam_index_t *idx, char *bamFileName, + char *position, bam_fetch_f callbackFunc, void *callbackData) +/* With the open bam file, return items the same way with the callbacks as with bamFetch() */ +/* except in this case use an already-open bam file and index (use bam_index_load and free() for */ +/* the index). It seems a little strange to pass the filename in with the open bam, but */ +/* it's just used to report errors. */ +{ +int chromId, start, end; +int ret = bam_parse_region(samfile->header, position, &chromId, &start, &end); +if (ret != 0 && startsWith("chr", position)) + ret = bam_parse_region(samfile->header, position+strlen("chr"), &chromId, &start, &end); +if (ret != 0) + // If the bam file does not cover the current chromosome, OK + return; +ret = bam_fetch(samfile->x.bam, idx, chromId, start, end, callbackData, callbackFunc); +if (ret != 0) + warn("bam_fetch(%s, %s (chromId=%d) failed (%d)", bamFileName, position, chromId, ret); +} + void bamFetch(char *fileOrUrl, char *position, bam_fetch_f callbackFunc, void *callbackData, samfile_t **pSamFile) /* Open the .bam file, fetch items in the seq:start-end position range, * and call callbackFunc on each bam item retrieved from the file plus callbackData. * This handles BAM files with "chr"-less sequence names, e.g. from Ensembl. * The pSamFile parameter is optional. If non-NULL it will be filled in, just for * the benefit of the callback function, with the open samFile. */ { char *bamFileName = NULL; samfile_t *fh = bamOpen(fileOrUrl, &bamFileName); boolean usingUrl = TRUE; usingUrl = (strstr(fileOrUrl, "tp://") || strstr(fileOrUrl, "https://")); if (pSamFile != NULL) *pSamFile = fh; -int chromId, start, end; -int ret = bam_parse_region(fh->header, position, &chromId, &start, &end); -if (ret != 0 && startsWith("chr", position)) - ret = bam_parse_region(fh->header, position+strlen("chr"), &chromId, &start, &end); -if (ret != 0) - // If the bam file does not cover the current chromosome, OK - return; #ifndef KNETFILE_HOOKS // Since samtools' url-handling code saves the .bai file to the current directory, // chdir to a trash directory before calling bam_index_load, then chdir back. char *runDir = getCurrentDir(); char *samDir = getSamDir(); if (usingUrl) setCurrentDir(samDir); #endif//ndef KNETFILE_HOOKS bam_index_t *idx = bam_index_load(bamFileName); #ifndef KNETFILE_HOOKS if (usingUrl) setCurrentDir(runDir); #endif//ndef KNETFILE_HOOKS if (idx == NULL) warn("bam_index_load(%s) failed.", bamFileName); else { - ret = bam_fetch(fh->x.bam, idx, chromId, start, end, callbackData, callbackFunc); - if (ret != 0) - warn("bam_fetch(%s, %s (chromId=%d) failed (%d)", bamFileName, position, chromId, ret); + bamFetchAlreadyOpen(fh, idx, bamFileName, position, callbackFunc, callbackData); free(idx); // Not freeMem, freez etc -- sam just uses malloc/calloc. } bamClose(&fh); } boolean bamIsRc(const bam1_t *bam) /* Return TRUE if alignment is on - strand. */ { const bam1_core_t *core = &bam->core; return (core->flag & BAM_FREVERSE); } void bamGetSoftClipping(const bam1_t *bam, int *retLow, int *retHigh, int *retClippedQLen) /* If retLow is non-NULL, set it to the number of "soft-clipped" (skipped) bases at * the beginning of the query sequence and quality; likewise for retHigh at end.