1 /*===========================================================================
2 *
3 * PUBLIC DOMAIN NOTICE
4 * National Center for Biotechnology Information
5 *
6 * This software/database is a "United States Government Work" under the
7 * terms of the United States Copyright Act. It was written as part of
8 * the author's official duties as a United States Government employee and
9 * thus cannot be copyrighted. This software/database is freely available
10 * to the public for use. The National Library of Medicine and the U.S.
11 * Government have not placed any restriction on its use or reproduction.
12 *
13 * Although all reasonable efforts have been taken to ensure the accuracy
14 * and reliability of the software and data, the NLM and the U.S.
15 * Government do not and cannot warrant the performance or results that
16 * may be obtained by using this software or data. The NLM and the U.S.
17 * Government disclaim all warranties, express or implied, including
18 * warranties of performance, merchantability or fitness for any particular
19 * purpose.
20 *
21 * Please cite the author in any work or product based on this material.
22 *
23 * ===========================================================================
24 *
25 */
26
27 #include <ncbi-vdb/NGS.hpp>
28 //#include <ngs-bam/ngs-bam.hpp>
29 #include <ngs/ErrorMsg.hpp>
30 #include <ngs/ReadCollection.hpp>
31 #include <ngs/Reference.hpp>
32 #include <ngs/AlignmentIterator.hpp>
33 #include <ngs/Alignment.hpp>
34
35 #include <math.h>
36 #include <iostream>
37
38 using namespace ngs;
39 using namespace std;
40
41 class AlignSliceTest
42 {
43 public:
44
run_common(ReadCollection & run,String refname,int start,int stop)45 static void run_common ( ReadCollection & run, String refname, int start, int stop )
46 {
47
48 String run_name = run.getName ();
49
50 // get requested reference
51 Reference ref = run.getReference ( refname );
52
53 // start iterator on requested range
54 long count = stop - start + 1;
55 AlignmentIterator it = ref.getAlignmentSlice ( start, count, Alignment::primaryAlignment );
56
57 long i;
58 for ( i = 0; it . nextAlignment (); ++ i )
59 {
60 cout << it.getReadId ()
61 << '\t' << it.getReferenceSpec ()
62 << '\t' << it.getAlignmentPosition ()
63 << '\t' << it.getLongCigar ( false ) // unclipped
64 << '\t' << it.getAlignedFragmentBases ()
65 << '\n';
66 }
67
68 cerr << "Read " << i << " alignments for " << run_name << '\n';
69 }
70
run_csra(String acc,String refname,int start,int stop)71 static void run_csra ( String acc, String refname, int start, int stop )
72 {
73 // open requested accession using SRA implementation of the API
74 ReadCollection run = ncbi::NGS::openReadCollection ( acc );
75 run_common ( run, refname, start, stop );
76 }
77
run_bam(String acc,String refname,int start,int stop)78 static void run_bam ( String acc, String refname, int start, int stop )
79 {
80 // open requested accession using example BAM implementation of the API
81 //ReadCollection run = NGS_BAM::openReadCollection ( acc );
82 //run_common ( run, refname, start, stop );
83 }
84
run(String acc,String refname,int start,int stop)85 static void run ( String acc, String refname, int start, int stop )
86 {
87 size_t dot = acc . find_last_of ( '.' );
88 if ( dot != string :: npos )
89 {
90 String extension = acc . substr ( dot );
91 if ( extension == ".bam" || extension == ".BAM" )
92 {
93 // run_bam ( acc, refname, start, stop );
94 return;
95 }
96 }
97
98 run_csra ( acc, refname, start, stop );
99 }
100
101
102 };
103
main(int argc,char const * argv[])104 int main ( int argc, char const *argv[] )
105 {
106 if ( argc != 5 )
107 {
108 cerr << "Usage: AlignSliceTest accession reference start stop\n";
109 }
110 else try
111 {
112 ncbi::NGS::setAppVersionString ( "AlignSliceTest.1.1.0" );
113 AlignSliceTest::run ( argv[1], argv[2], atoi ( argv[3] ), atoi ( argv[4] ) );
114 return 0;
115 }
116 catch ( ErrorMsg & x )
117 {
118 cerr << x.toString () << '\n';
119 }
120 catch ( exception & x )
121 {
122 cerr << x.what () << '\n';
123 }
124 catch ( ... )
125 {
126 cerr << "unknown exception\n";
127 }
128
129 return 10;
130 }
131