incbwt/read_bwt.cpp

   1 #include <cstdlib>
   2 #include <fstream>
   3 #include <iostream>
   4
   5 #include "rlcsa.h"
   6
   7
   8 using namespace CSA;
   9
  10
  11 /*
  12   This program writes run-length encoded PLCP of the collection into a file.
  13 */
  14
  15
  16 usint
  17 countRuns(usint prev, uchar* buffer, usint length)
  18 {
  19   usint runs = 0;
  20
  21   for(usint i = 0; i < length; i++)
  22   {
  23     if(buffer[i] != prev)
  24     {
  25       prev = buffer[i];
  26       if(buffer[i] != 0) { runs++; }
  27     }
  28   }
  29
  30   return runs;
  31 }
  32
  33
  34 int
  35 main(int argc, char** argv)
  36 {
  37   std::cout << "RLCSA to BWT converter" << std::endl;
  38   if(argc < 2)
  39   {
  40     std::cout << "Usage: read_bwt base_name [buffer_size]" << std::endl;
  41     return 1;
  42   }
  43
  44   std::string base_name = argv[1];
  45   std::string bwt_name = base_name + ".bwt";
  46   std::cout << "BWT: " << bwt_name << std::endl;
  47   std::ofstream bwt_file(bwt_name.c_str(), std::ios_base::binary);
  48   if(!bwt_file)
  49   {
  50     std::cerr << "Error creating BWT file!" << std::endl;
  51     return 2;
  52   }
  53   std::cout << std::endl;
  54
  55   RLCSA rlcsa(base_name);
  56   clock_t start = clock();
  57   usint buffer_size = 0;
  58   if(argc > 2) { buffer_size = atoi(argv[2]); }
  59   usint n = rlcsa.getSize() + rlcsa.getNumberOfSequences();
  60
  61   usint runs = 0, prev = CHARS;
  62   if(buffer_size > 0)
  63   {
  64     for(usint i = 0; i < n; i += buffer_size)
  65     {
  66       pair_type range(i, std::min(i + buffer_size - 1, n - 1));
  67       uchar* bwt = rlcsa.readBWT(range);
  68       if(bwt != 0)
  69       {
  70         runs += countRuns(prev, bwt, length(range));
  71         prev = bwt[length(range) - 1];
  72         bwt_file.write((char*)bwt, length(range));
  73         delete[] bwt;
  74       }
  75     }
  76   }
  77   else
  78   {
  79     uchar* bwt = rlcsa.readBWT();
  80     if(bwt != 0)
  81     {
  82       runs = countRuns(prev, bwt, n);
  83       bwt_file.write((char*)bwt, n);
  84       delete[] bwt;
  85     }
  86   }
  87
  88   clock_t stop = clock();
  89   double time = ((stop - start) / (double)CLOCKS_PER_SEC);
  90   double megabytes = n / (double)MEGABYTE;
  91   std::cout << megabytes << " megabytes in " << time << " seconds (" << (megabytes / time) << " MB/s)" << std::endl;
  92   std::cout << std::endl;
  93
  94   // Testing direct reporting of the number of runs.
  95   // This is as expensive as reading the BWT.
  96   std::cout << "Number of runs: " << runs << std::endl;
  97   runs = rlcsa.countRuns();
  98   std::cout << "Number of runs (direct count): " << runs << std::endl;
  99   std::cout << std::endl;
 100
 101   bwt_file.close();
 102   return 0;
 103 }