Skip to content

Commit

Permalink
now prints alignments to contigs which do not appear in the vcf as th…
Browse files Browse the repository at this point in the history
…ey were (previously the RNAME is empty for such contigs); also prints the contigs in the header of lifted sam
  • Loading branch information
milkschen committed Feb 25, 2020
1 parent 30c3bb1 commit eea21c7
Show file tree
Hide file tree
Showing 2 changed files with 11 additions and 2 deletions.
8 changes: 7 additions & 1 deletion liftover.cpp
Original file line number Diff line number Diff line change
Expand Up @@ -143,6 +143,12 @@ void lift_run(lift_opts args) {
for (auto i = 0; i < contig_names.size(); ++i) {
fprintf(out_sam_fp, "@SQ\tSN:%s\tLN:%ld\n", contig_names[i].data(), contig_reflens[i]);
}
for (auto i = 0; i < hdr->n_targets; i++){
// if a contig is not found in vcf/lft, print it as it was before liftover
if (std::find(contig_names.begin(), contig_names.end(), hdr->target_name[i]) == contig_names.end()){
fprintf(out_sam_fp, "@SQ\tSN:%s\tLN:%u\n", hdr->target_name[i], hdr->target_len[i]);
}
}
/* This is only supported in the latest dev release of htslib as of July 22.
* add back in at next htslib release
kstring_t s = KS_INITIALIZE;
Expand Down Expand Up @@ -278,7 +284,7 @@ int main(int argc, char** argv) {
exit(1);
}
if (args.haplotype != "0" && args.haplotype != "1"){
fprintf(stderr, "invalid haplotype %s\n", args.haplotype);
fprintf(stderr, "invalid haplotype %s\n", args.haplotype.c_str());
exit(1);
}
if (!strcmp(argv[optind], "lift")) {
Expand Down
5 changes: 4 additions & 1 deletion liftover.hpp
Original file line number Diff line number Diff line change
Expand Up @@ -368,7 +368,10 @@ class LiftMap {
std::string get_other_name(std::string n) {
if (name_map.find(n) != name_map.end()) {
return name_map[n];
} else return "";
}
// return original name if mapping is not available
else
return n;
}

// converts CIGAR string of s2 alignment to corresponding CIGAR string for the s1 alignment
Expand Down

0 comments on commit eea21c7

Please sign in to comment.