C++使用htslib庫讀入和寫出bam文件的實例
有時候我們需要使用C++處理bam文件,比如取出read1或者read2等符合特定條件的序列,根據(jù)cigar值對序列指定位置的堿基進行統(tǒng)計或者對序列進行處理并輸出等,這時我們可以使用htslib庫。htslib可以用來處理SAM, BAM,CRAM 和VCF文件,是samtools、bcftools的核心庫。
#include <stdio.h> #include <stdlib.h> #include <htslib/sam.h> using namespace std; #define bam_is_read1(b) (((b)->core.flag&BAM_FREAD1) != 0) uint8_t Base[16] = {0,65,67,0,71,0,0,0,84,0,0,0,0,0,0,78}; int main(int argc, char **argv) { bam_hdr_t *header; bam1_t *aln = bam_init1(); samFile *in = sam_open(argv[1], "r"); htsFile *outR1 = hts_open(argv[2], "wb"); header = sam_hdr_read(in); if (sam_hdr_write(outR1, header) < 0) { fprintf(stderr, "Error writing output.\n"); exit(-1); } uint8_t *seq; int32_t lseq; uint32_t *cigar; char* qname; while (sam_read1(in, header, aln) >= 0) { if (bam_is_read1(aln)){ sam_write1(outR1, header, aln); } else { seq = bam_get_seq(aln); lseq = aln->core.l_qseq; qname = bam_get_qname(aln); printf("%s\n",qname); cigar = bam_get_cigar(aln); for(int i=0; i < aln->core.n_cigar;++i){ int icigar = cigar[i]; printf("%d%d\n",bam_cigar_op(icigar),bam_cigar_oplen(icigar)); } for(int i=0; i < lseq;++i){ printf("%c", Base[bam_seqi(seq, i)]); } printf("\n"); } } sam_close(in); sam_close(outR1); }
cigar值存儲形式
32位int,通過bam_get_cigar獲得地址,aln->core.n_cigar返回cigar operation的個數(shù)
•低 4位存儲cigar operation;通過函數(shù)bam_cigar_op()獲得operation
•高28位存儲cigar值的長度;通過函數(shù),bam_cigar_oplen獲得
seq存儲形式
8位int,4位存儲一個堿基,1,2,4,8,15分別代表A、C、G、T、N,高四位存儲坐標(biāo)數(shù)較小的堿基,可通過bam_seqi(seq,i)遍歷。
以上這篇C++使用htslib庫讀入和寫出bam文件的實例就是小編分享給大家的全部內(nèi)容了,希望能給大家一個參考,也希望大家多多支持我們。
上一篇:基于C++ bitset常用函數(shù)及運算符(詳解)
欄 目:C語言
本文標(biāo)題:C++使用htslib庫讀入和寫出bam文件的實例
本文地址:http://mengdiqiu.com.cn/a1/Cyuyan/1040.html
您可能感興趣的文章
- 04-02func函數(shù)+在C語言 func函數(shù)在c語言中
- 04-02c語言沒有round函數(shù) round c語言
- 01-10使用OpenGL實現(xiàn)3D立體顯示的程序代碼
- 01-10深入理解C++中常見的關(guān)鍵字含義
- 01-10使用C++實現(xiàn)全排列算法的方法詳解
- 01-10深入Main函數(shù)中的參數(shù)argc,argv的使用詳解
- 01-10c++中inline的用法分析
- 01-10用C++實現(xiàn)DBSCAN聚類算法
- 01-10全排列算法的非遞歸實現(xiàn)與遞歸實現(xiàn)的方法(C++)
- 01-10C++大數(shù)模板(推薦)


閱讀排行
本欄相關(guān)
- 04-02c語言函數(shù)調(diào)用后清空內(nèi)存 c語言調(diào)用
- 04-02func函數(shù)+在C語言 func函數(shù)在c語言中
- 04-02c語言的正則匹配函數(shù) c語言正則表達
- 04-02c語言用函數(shù)寫分段 用c語言表示分段
- 04-02c語言中對數(shù)函數(shù)的表達式 c語言中對
- 04-02c語言編寫函數(shù)冒泡排序 c語言冒泡排
- 04-02c語言沒有round函數(shù) round c語言
- 04-02c語言分段函數(shù)怎么求 用c語言求分段
- 04-02C語言中怎么打出三角函數(shù) c語言中怎
- 04-02c語言調(diào)用函數(shù)求fibo C語言調(diào)用函數(shù)求
隨機閱讀
- 08-05DEDE織夢data目錄下的sessions文件夾有什
- 08-05織夢dedecms什么時候用欄目交叉功能?
- 01-10delphi制作wav文件的方法
- 01-10SublimeText編譯C開發(fā)環(huán)境設(shè)置
- 04-02jquery與jsp,用jquery
- 01-10使用C語言求解撲克牌的順子及n個骰子
- 08-05dedecms(織夢)副欄目數(shù)量限制代碼修改
- 01-11Mac OSX 打開原生自帶讀寫NTFS功能(圖文
- 01-11ajax實現(xiàn)頁面的局部加載
- 01-10C#中split用法實例總結(jié)