Skip to content

Commit

Permalink
Also accept SAM headers from <baseprefix>.dict if present
Browse files Browse the repository at this point in the history
Now outputs SAM headers from: -H options; .hdr file stored with the
index files if present, or otherwise a .dict file if that is present;
the basic hard-coded `@HD`/`@SQ` default headers.
  • Loading branch information
jmarshall committed Mar 18, 2024
1 parent 46dee88 commit d55a779
Show file tree
Hide file tree
Showing 2 changed files with 16 additions and 4 deletions.
6 changes: 4 additions & 2 deletions bwa.1
Original file line number Diff line number Diff line change
Expand Up @@ -85,7 +85,7 @@ Identifies reference sequences that are alternative alleles, and provides
alignment information that can be used to lift ALT hits over to corresponding
positions within the primary assembly.
.TP
.IR prefix .hdr
.IR prefix ".hdr or " baseprefix .dict
Provides SAM headers to be used as the mapping output's headers, instead of
the default simple `@HD' and `@SQ' headers. This can be used to add descriptive
fields to the `@SQ' headers output: species, assembly, alternative names (which
Expand Down Expand Up @@ -325,7 +325,9 @@ added via
.B -H
include any `@SQ' headers, the
.IR db.prefix .hdr
file (if it exists) will be ignored. [null]
and
.IR db.baseprefix .dict
files (if either exist) will be ignored. [null]
.TP
.BI -o \ FILE
Write the output SAM file to
Expand Down
14 changes: 12 additions & 2 deletions bwa.c
Original file line number Diff line number Diff line change
Expand Up @@ -484,12 +484,22 @@ void bwa_print_sam_hdr2(const bntseq_t *bns, const char *prefix, const char *hdr
kstring_t str = { 0, 0, NULL };
char *bns_hdr = NULL;

// Ignore index .hdr file entirely if the user's headers provide @SQ lines
// Ignore index .hdr/.dict file entirely if the user's headers provide @SQ lines
if (!has_SQ(hdr_line)) {
FILE *fp;
// Otherwise read the .hdr file if present
ksprintf(&str, "%s.hdr", prefix);
if ((fp = fopen(str.s, "r")) != NULL) {
fp = fopen(str.s, "r");
if (!fp) { // Otherwise also try .dict if present
size_t l;
str.l -= 4;
if (str.l >= 3 && strncmp(&str.s[str.l-3], ".gz", 3) == 0) str.l -= 3;
for (l = str.l; l > 0 && str.s[l-1] != '/'; l--)
if (str.s[l-1] == '.') { str.l = l-1; break; }
ksprintf(&str, ".dict");
fp = fopen(str.s, "r");
}
if (fp) {
int c, n_SQ;
str.l = 0;
while ((c = getc(fp)) != EOF)
Expand Down

0 comments on commit d55a779

Please sign in to comment.