X-Git-Url: https://git.donarmstrong.com/?a=blobdiff_plain;f=bam_import.c;h=29516c92ff42593313caf2833c38213d62ebdc52;hb=80658e4d152b53bd55927c432ceece7702ab49d7;hp=2b4fb7e03f956e6ea35a2216c5c853c414fa0aa4;hpb=ed9cedc66455cf33331261e0d230c508411b3df6;p=samtools.git diff --git a/bam_import.c b/bam_import.c index 2b4fb7e..29516c9 100644 --- a/bam_import.c +++ b/bam_import.c @@ -116,7 +116,7 @@ static bam_header_t *hash2header(const kh_ref_t *hash) bam_header_t *sam_header_read2(const char *fn) { bam_header_t *header; - int c, dret, ret; + int c, dret, ret, error = 0; gzFile fp; kstream_t *ks; kstring_t *str; @@ -135,6 +135,10 @@ bam_header_t *sam_header_read2(const char *fn) ks_getuntil(ks, 0, str, &dret); len = atoi(str->s); k = kh_put(ref, hash, s, &ret); + if (ret == 0) { + fprintf(stderr, "[sam_header_read2] duplicated sequence name: %s\n", s); + error = 1; + } kh_value(hash, k) = (uint64_t)len<<32 | i; if (dret != '\n') while ((c = ks_getc(ks)) != '\n' && c != -1); @@ -143,6 +147,7 @@ bam_header_t *sam_header_read2(const char *fn) gzclose(fp); free(str->s); free(str); fprintf(stderr, "[sam_header_read2] %d sequences loaded.\n", kh_size(hash)); + if (error) return 0; header = hash2header(hash); kh_destroy(ref, hash); return header; @@ -163,21 +168,29 @@ static inline void parse_error(int64_t n_lines, const char * __restrict msg) } static inline void append_text(bam_header_t *header, kstring_t *str) { - int x = header->l_text, y = header->l_text + str->l + 2; // 2 = 1 byte dret + 1 byte null + size_t x = header->l_text, y = header->l_text + str->l + 2; // 2 = 1 byte dret + 1 byte null kroundup32(x); kroundup32(y); - if (x < y) header->text = (char*)realloc(header->text, y); + if (x < y) + { + header->n_text = y; + header->text = (char*)realloc(header->text, y); + if ( !header->text ) + { + fprintf(stderr,"realloc failed to alloc %ld bytes\n", y); + abort(); + } + } + // Sanity check + if ( header->l_text+str->l+1 >= header->n_text ) + { + fprintf(stderr,"append_text FIXME: %ld>=%ld, x=%ld,y=%ld\n", header->l_text+str->l+1,header->n_text,x,y); + abort(); + } strncpy(header->text + header->l_text, str->s, str->l+1); // we cannot use strcpy() here. header->l_text += str->l + 1; header->text[header->l_text] = 0; } -int sam_header_parse_rg(bam_header_t *h) -{ - if (h->dict == 0) h->dict = sam_header_parse2(h->text); - if (h->rg2lib) h->rg2lib = sam_header2tbl(h->dict, "RG", "ID", "LB"); - return sam_tbl_size(h->rg2lib); -} - int sam_header_parse(bam_header_t *h) { char **tmp; @@ -186,10 +199,14 @@ int sam_header_parse(bam_header_t *h) h->n_targets = 0; h->target_len = 0; h->target_name = 0; if (h->l_text < 3) return 0; if (h->dict == 0) h->dict = sam_header_parse2(h->text); - h->target_name = sam_header2list(h->dict, "SQ", "SN", &h->n_targets); + tmp = sam_header2list(h->dict, "SQ", "SN", &h->n_targets); if (h->n_targets == 0) return 0; + h->target_name = calloc(h->n_targets, sizeof(void*)); + for (i = 0; i < h->n_targets; ++i) + h->target_name[i] = strdup(tmp[i]); + free(tmp); tmp = sam_header2list(h->dict, "SQ", "LN", &h->n_targets); - h->target_len = (uint32_t*)calloc(h->n_targets, 4); + h->target_len = calloc(h->n_targets, 4); for (i = 0; i < h->n_targets; ++i) h->target_len[i] = atoi(tmp[i]); free(tmp); @@ -410,6 +427,27 @@ int sam_read1(tamFile fp, bam_header_t *header, bam1_t *b) memcpy(s, str->s + 5, str->l - 5); s[str->l - 5] = 0; doff += size; + } else if (type == 'B') { + int32_t n = 0, Bsize, k = 0, size; + char *p; + if (str->l < 8) parse_error(fp->n_lines, "too few values in aux type B"); + Bsize = bam_aux_type2size(str->s[5]); // the size of each element + for (p = (char*)str->s + 6; *p; ++p) // count the number of elements in the array + if (*p == ',') ++n; + p = str->s + 7; // now p points to the first number in the array + size = 6 + Bsize * n; // total number of bytes allocated to this tag + s = alloc_data(b, doff + 6 * Bsize * n) + doff; // allocate memory + *s++ = 'B'; *s++ = str->s[5]; + memcpy(s, &n, 4); s += 4; // write the number of elements + if (str->s[5] == 'c') while (p < str->s + str->l) ((int8_t*)s)[k++] = (int8_t)strtol(p, &p, 0), ++p; + else if (str->s[5] == 'C') while (p < str->s + str->l) ((uint8_t*)s)[k++] = (uint8_t)strtol(p, &p, 0), ++p; + else if (str->s[5] == 's') while (p < str->s + str->l) ((int16_t*)s)[k++] = (int16_t)strtol(p, &p, 0), ++p; // FIXME: avoid unaligned memory + else if (str->s[5] == 'S') while (p < str->s + str->l) ((uint16_t*)s)[k++] = (uint16_t)strtol(p, &p, 0), ++p; + else if (str->s[5] == 'i') while (p < str->s + str->l) ((int32_t*)s)[k++] = (int32_t)strtol(p, &p, 0), ++p; + else if (str->s[5] == 'I') while (p < str->s + str->l) ((uint32_t*)s)[k++] = (uint32_t)strtol(p, &p, 0), ++p; + else if (str->s[5] == 'f') while (p < str->s + str->l) ((float*)s)[k++] = (float)strtod(p, &p), ++p; + else parse_error(fp->n_lines, "unrecognized array type"); + s += Bsize * n; doff += size; } else parse_error(fp->n_lines, "unrecognized type"); if (dret == '\n' || dret == '\r') break; }