test/bam_pair_buffer.cc

Code
Comments
Other
Rev Date Author Line
4057 31 Mar 21 peter 1 // $Id$
4057 31 Mar 21 peter 2 //
4057 31 Mar 21 peter 3 // Copyright (C) 2021 Peter Johansson
4057 31 Mar 21 peter 4 //
4057 31 Mar 21 peter 5 // This program is free software; you can redistribute it and/or modify
4057 31 Mar 21 peter 6 // it under the terms of the GNU General Public License as published by
4057 31 Mar 21 peter 7 // the Free Software Foundation; either version 3 of the License, or
4057 31 Mar 21 peter 8 // (at your option) any later version.
4057 31 Mar 21 peter 9 //
4057 31 Mar 21 peter 10 // This program is distributed in the hope that it will be useful, but
4057 31 Mar 21 peter 11 // WITHOUT ANY WARRANTY; without even the implied warranty of
4057 31 Mar 21 peter 12 // MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
4057 31 Mar 21 peter 13 // General Public License for more details.
4057 31 Mar 21 peter 14 //
4057 31 Mar 21 peter 15 // You should have received a copy of the GNU General Public License
4057 31 Mar 21 peter 16 // along with this program. If not, see <http://www.gnu.org/licenses/>.
4057 31 Mar 21 peter 17
4057 31 Mar 21 peter 18 #include <config.h>
4057 31 Mar 21 peter 19
4057 31 Mar 21 peter 20 #include "Suite.h"
4057 31 Mar 21 peter 21
4057 31 Mar 21 peter 22 #if YAT_HAVE_HTSLIB
4057 31 Mar 21 peter 23 #include "yat/omic/algorithm.h"
4057 31 Mar 21 peter 24 #include "yat/omic/BamFile.h"
4057 31 Mar 21 peter 25 #include "yat/omic/BamPairBuffer.h"
4057 31 Mar 21 peter 26 #include "yat/omic/BamReadFilter.h"
4057 31 Mar 21 peter 27 #include "yat/omic/BamPairIterator2.h"
4057 31 Mar 21 peter 28 #include "yat/omic/BamReadIterator.h"
4057 31 Mar 21 peter 29 #endif
4057 31 Mar 21 peter 30
4057 31 Mar 21 peter 31 #include <boost/iterator/filter_iterator.hpp>
4057 31 Mar 21 peter 32
4057 31 Mar 21 peter 33 #include <iostream>
4057 31 Mar 21 peter 34
4057 31 Mar 21 peter 35 using namespace theplu::yat;
4057 31 Mar 21 peter 36
4057 31 Mar 21 peter 37 void test1(test::Suite& suite);
4057 31 Mar 21 peter 38
4057 31 Mar 21 peter 39 int main(int argc, char* argv[])
4057 31 Mar 21 peter 40 {
4057 31 Mar 21 peter 41   test::Suite suite(argc, argv, true);
4057 31 Mar 21 peter 42
4057 31 Mar 21 peter 43 #ifdef YAT_HAVE_HTSLIB
4057 31 Mar 21 peter 44   test1(suite);
4057 31 Mar 21 peter 45 #endif
4057 31 Mar 21 peter 46   return suite.return_value();
4057 31 Mar 21 peter 47 }
4057 31 Mar 21 peter 48
4057 31 Mar 21 peter 49
4057 31 Mar 21 peter 50 #ifdef YAT_HAVE_HTSLIB
4057 31 Mar 21 peter 51 using namespace omic;
4057 31 Mar 21 peter 52
4057 31 Mar 21 peter 53
4057 31 Mar 21 peter 54 template<typename InputIterator, typename OutputIterator>
4057 31 Mar 21 peter 55 void process(InputIterator first, InputIterator last,
4057 31 Mar 21 peter 56              OutputIterator out)
4057 31 Mar 21 peter 57 {
4057 31 Mar 21 peter 58   BamPairBuffer<OutputIterator> buffer(out);
4057 31 Mar 21 peter 59   BamPairIterator2<InputIterator> end(last, last);
4057 31 Mar 21 peter 60   for (BamPairIterator2<InputIterator> iter(first, last); iter!=end; ++iter)
4057 31 Mar 21 peter 61     buffer.push(iter->first(), iter->second());
4057 31 Mar 21 peter 62 }
4057 31 Mar 21 peter 63
4057 31 Mar 21 peter 64
4057 31 Mar 21 peter 65 void test1(test::Suite& suite)
4057 31 Mar 21 peter 66 {
4057 31 Mar 21 peter 67   std::string file = "../../data/foo.sorted.bam";
4057 31 Mar 21 peter 68
4057 31 Mar 21 peter 69   InBamFile is(file);
4057 31 Mar 21 peter 70   BamReadIterator iter(is);
4057 31 Mar 21 peter 71   BamReadIterator end;
4057 31 Mar 21 peter 72   std::vector<BamRead> result;
4057 31 Mar 21 peter 73   BamReadFilter filter(0, 0, 2);
4057 31 Mar 21 peter 74   process(boost::make_filter_iterator(filter, iter, end),
4057 31 Mar 21 peter 75           boost::make_filter_iterator(filter, end),
4057 31 Mar 21 peter 76           std::back_inserter(result));
4057 31 Mar 21 peter 77   is.close();
4057 31 Mar 21 peter 78
4057 31 Mar 21 peter 79   if (result.empty()) {
4057 31 Mar 21 peter 80     suite.add(false);
4057 31 Mar 21 peter 81     suite.err() << "error: empty result\n";
4057 31 Mar 21 peter 82   }
4057 31 Mar 21 peter 83
4057 31 Mar 21 peter 84   for (size_t i=0; i<result.size(); ++i) {
4057 31 Mar 21 peter 85     if ((result[i].flag() & 2) == 0) {
4057 31 Mar 21 peter 86       suite.add(false);
4057 31 Mar 21 peter 87       suite.err() << "error: bit x02 not set in flag\n";
4057 31 Mar 21 peter 88     }
4057 31 Mar 21 peter 89   }
4057 31 Mar 21 peter 90
4057 31 Mar 21 peter 91   if (!std::is_sorted(result.begin(), result.end(), BamLessPos())) {
4057 31 Mar 21 peter 92     suite.add(false);
4057 31 Mar 21 peter 93     suite.err() << "error: result is not sorted\n";
4057 31 Mar 21 peter 94     BamLessPos less;
4057 31 Mar 21 peter 95     for (size_t i=0; i+1<result.size(); ++i)
4057 31 Mar 21 peter 96       if (less(result[i+1], result[i]))
4057 31 Mar 21 peter 97         suite.err() << "position of " << result[i].name() << ", "
4057 31 Mar 21 peter 98                     << result[i].pos() << ", is greater than "
4057 31 Mar 21 peter 99                     << "position of " << result[i+1].name() << ", "
4057 31 Mar 21 peter 100                     << result[i+1].pos() << ".\n";
4057 31 Mar 21 peter 101   }
4057 31 Mar 21 peter 102
4057 31 Mar 21 peter 103 }
4057 31 Mar 21 peter 104 #endif

Generated on Sat Oct 14 02:20:41 2023 (UTC) by svndigest 0.11.2pre (dev build 1672)