dcover/distribute.hpp

   1 /*
   2  * distribute.hpp
   3  *
   4  * Copyright 2004 Juha K"arkk"ainen <juha.karkkainen@cs.helsinki.fi>
   5  *
   6  */
   7
   8 #ifndef DISTRIBUTE_HPP
   9 #define DISTRIBUTE_HPP
  10
  11 #include "dcsuffixsort.hpp"
  12
  13 #include <iterator>
  14 #include <algorithm>
  15 #include <utility>
  16
  17
  18 //----------------------------------------------------------------------
  19 // suffix_lcp
  20 //
  21 // Compute the length of the longest common prefix between two
  22 // suffixes a and b.
  23 // str_end = end of the text
  24 // minlcp = length of the already known common prefix
  25 // maxlcp = stop when the length reaches maxlcp
  26 //----------------------------------------------------------------------
  27 template <class StringIterator, class LcpType>
  28 LcpType
  29 suffix_lcp(StringIterator a, StringIterator b, StringIterator str_end,
  30            LcpType minlcp, LcpType maxlcp)
  31 {
  32   if (a < b) swap(a,b);   // make a the shorter suffix
  33   if ((str_end-a) > maxlcp) str_end = a + maxlcp;
  34   std::pair<StringIterator,StringIterator> mm =
  35               std::mismatch(a+minlcp, str_end, b+minlcp);
  36   //std::cout << " lcp=" << mm.first-a << "\n";
  37   return static_cast<LcpType>(mm.first-a);
  38 }
  39
  40
  41 //----------------------------------------------------------------------
  42 // class suffix_distributor
  43 //----------------------------------------------------------------------
  44 template <class StringIterator, class SuffixIterator,
  45           class DCSampler, class RankIterator>
  46 class suffix_distributor
  47 {
  48 private:
  49   //typedef typename std::iterator_traits<StringIterator>::difference_type
  50   //  lcptype;
  51   typedef short lcptype;
  52   typedef typename
  53     std::iterator_traits<SuffixIterator>::difference_type difftype;
  54
  55
  56   StringIterator str, str_end;
  57   SuffixIterator pivots, pivots_end;
  58   std::vector<lcptype> lcps;
  59   const DCSampler& dcsampler;
  60   RankIterator dcranks;
  61   lcptype maxlcp;
  62
  63   typedef typename std::vector<lcptype>::iterator lcpiterator;
  64
  65 public:
  66
  67   suffix_distributor(StringIterator str_, StringIterator str_end_,
  68                      SuffixIterator pivots_, SuffixIterator pivots_end_,
  69                      const DCSampler& dcsampler_, RankIterator dcranks_)
  70     : str(str_), str_end(str_end_), pivots(pivots_), pivots_end(pivots_end_),
  71       lcps(pivots_end-pivots), dcsampler(dcsampler_), dcranks(dcranks_),
  72     maxlcp(dcsampler.period())
  73   {
  74     compute_lcps();
  75   }
  76
  77   void reset(SuffixIterator pivots_, SuffixIterator pivots_end_)
  78   {
  79     pivots=pivots_;
  80     pivots_end=pivots_end_;
  81     lcps.resize(pivots_end-pivots);
  82
  83     compute_lcps();
  84   }
  85
  86
  87   difftype nbuckets() const{ return pivots_end-pivots+1; }
  88
  89   difftype find_bucket(StringIterator suf) const {
  90     return binary_search(suf);
  91     //return i - pivots;
  92   }
  93
  94   void
  95   print() const {
  96 #if DEBUG > 1
  97     std::cout << "distributor memory="
  98               << lcps.size()*sizeof(lcptype)
  99               << std::endl;
 100 #if DEBUG > 2
 101     std::cout << "pivots and lcps:\n";
 102     for (SuffixIterator i = pivots; i != pivots_end; ++i) {
 103       std::cout << std::setw(4) << i-pivots
 104                 << std::setw(4) << int(lcps[i-pivots])
 105                 << std::setw(9) << *i << " ";
 106       std::copy(str+*i, std::min(str_end, str+*i+40),
 107                 std::ostream_iterator<char>(std::cout,""));
 108       std::cout << "\n";
 109     }
 110     std::cout << "distributor end\n";
 111 #endif
 112 #endif
 113   }
 114
 115
 116 private:
 117
 118   bool dc_less_or_equal(difftype a, difftype b) const {
 119     std::pair<difftype,difftype> rep = dcsampler.pair(a, b);
 120     return dcranks[rep.first] <= dcranks[rep.second];
 121   }
 122
 123   //----------------------------------------------------------------------
 124   // compute_lcps_bounded
 125   //
 126   // Compute all the lcp values in a subrange bounded by
 127   // actual pivot suffixes.
 128   //----------------------------------------------------------------------
 129   void
 130   compute_lcps_bounded(SuffixIterator left_suf, SuffixIterator right_suf,
 131                       lcpiterator left_lcp,
 132                       lcptype minlcp, lcptype maxlcp)
 133   {
 134     if (left_suf == right_suf) return;
 135     difftype midpoint = (right_suf-left_suf)/2;
 136     //std::cout << "  midpoint of [0," << right_suf-left_suf
 137     //          << ")=" << midpoint << "\n";
 138
 139     SuffixIterator mid_suf = left_suf + midpoint;
 140     lcpiterator mid_lcp = left_lcp + midpoint;
 141
 142     lcptype llcp_left = suffix_lcp(str+*(left_suf-1), str+*mid_suf,
 143                                    str_end, minlcp, maxlcp);
 144     lcptype llcp_right = suffix_lcp(str+*mid_suf, str+*right_suf,
 145                                     str_end, minlcp, maxlcp);
 146     *mid_lcp = llcp_right - llcp_left;
 147     //std::cout << *mid_lcp << std::endl;
 148
 149     compute_lcps_bounded(left_suf, mid_suf, left_lcp, llcp_left, maxlcp);
 150     compute_lcps_bounded(mid_suf+1, right_suf, mid_lcp+1, llcp_right, maxlcp);
 151   }
 152
 153
 154   //----------------------------------------------------------------------
 155   // compute_lcps
 156   //
 157   // Compute all the lcp values in the full range.
 158   //----------------------------------------------------------------------
 159   void
 160   compute_lcps()
 161   {
 162     if (pivots == pivots_end) return;
 163     difftype centerpoint = (pivots_end - pivots)/2;
 164     //std::cout << "midpoint of [0," << (pivots_end - pivots) << ")="
 165     //      << centerpoint << "\n";
 166     lcps[centerpoint] = 0;
 167
 168     // unbounded ranges on the left
 169     difftype rightpoint = centerpoint;
 170     while (rightpoint > 0) {
 171       difftype midpoint = rightpoint/2;
 172       //std::cout << "midpoint of [0," << rightpoint << ")=" << midpoint << "\n";
 173       lcptype llcp = suffix_lcp(str+pivots[midpoint],
 174                                str+pivots[rightpoint],
 175                                str_end, static_cast<lcptype>(0), maxlcp);
 176       lcps[midpoint] = llcp;
 177       //std::cout << llcp << std::endl;
 178       compute_lcps_bounded(pivots+midpoint+1, pivots+rightpoint,
 179                            lcps.begin()+midpoint+1, llcp, maxlcp);
 180       rightpoint = midpoint;
 181     }
 182
 183     // unbounded ranges on the right
 184     SuffixIterator left_suf = pivots + centerpoint + 1;
 185     lcpiterator left_lcp = lcps.begin() + centerpoint + 1;
 186     while (left_suf != pivots_end) {
 187       difftype midpoint = (pivots_end - left_suf)/2;
 188       //std::cout << "midpoint of [" << (left_suf-pivots) << "," <<
 189       //  (pivots_end-pivots) << ")=" << (left_suf-pivots)+midpoint << "\n";
 190       lcptype llcp = suffix_lcp(str+left_suf[-1], str+left_suf[midpoint],
 191                                str_end, static_cast<lcptype>(0), maxlcp);
 192       left_lcp[midpoint] = -llcp;
 193       //std::cout << -llcp << std::endl;
 194       compute_lcps_bounded(left_suf, left_suf+midpoint, left_lcp,
 195                            llcp, maxlcp);
 196       left_suf += midpoint + 1;
 197       left_lcp += midpoint + 1;
 198     }
 199
 200   }
 201
 202
 203
 204 #ifndef BINSEARCH
 205
 206 //----------------------------------------------------------------------
 207 // binary_search
 208 //----------------------------------------------------------------------
 209 difftype
 210 //SuffixIterator
 211 binary_search(StringIterator query) const
 212 {
 213   //SuffixIterator left = pivots;
 214   //SuffixIterator right = pivots_end;
 215   //SuffixIterator mid;
 216
 217   difftype left = 0;
 218   difftype right = pivots_end-pivots;
 219   difftype mid;
 220
 221   lcptype minlcp = 0;
 222   lcptype lcp_diff = 0;
 223   lcptype lcp_left, lcp_right;
 224
 225  in_the_middle:
 226   while (true) {
 227     /*
 228     std::cout << "middle minlcp=" << int(minlcp)
 229               << " lcp_diff=" << int(lcp_diff)
 230               << " mid=" << left + (right-left)/2
 231               << "\n";
 232     */
 233     if (left == right) return left;
 234     mid = left + (right-left)/2;
 235     if (str+pivots[mid] == query) return mid+1;
 236     lcptype lcp_mid = lcps[mid];
 237     if (lcp_mid == 0) {
 238       lcptype newlcp = suffix_lcp(query, str+pivots[mid], str_end,
 239                                   minlcp, maxlcp);
 240       if (newlcp == maxlcp) {
 241         lcp_left = lcp_right = minlcp;
 242         goto fullmatch;
 243       }
 244       if (query+newlcp == str_end ||
 245           query[newlcp] < str[pivots[mid]+newlcp]) {
 246         right = mid;
 247         lcp_diff = newlcp - minlcp;
 248         if (lcp_diff != 0) goto closer_to_right;
 249       } else {
 250         left = mid+1;
 251         lcp_diff = minlcp - newlcp;
 252         if (lcp_diff != 0) goto closer_to_left;
 253       }
 254     } else {
 255       if (lcp_mid < lcp_diff) {
 256         left = mid+1;
 257       } else {
 258         right = mid;
 259       }
 260     }
 261   }
 262
 263  closer_to_left:
 264   while (true) {
 265     /*
 266     std::cout << "left minlcp=" << int(minlcp)
 267               << " lcp_diff=" << int(lcp_diff)
 268               << " mid=" << left + (right-left)/2
 269               << "\n";
 270     */
 271     if (left == right) return left;
 272     mid = left + (right-left)/2;
 273     if (str+pivots[mid] == query) return mid+1;
 274     lcptype lcp_mid = lcps[mid];
 275     if (lcp_mid > lcp_diff) {
 276       if (lcp_mid < 0) {
 277         lcp_diff -= lcp_mid;
 278         minlcp -= lcp_mid;
 279       }
 280       right = mid;
 281     } else {
 282       if (lcp_mid == lcp_diff) {
 283         lcptype lcp = minlcp - lcp_diff;
 284         lcptype newlcp = suffix_lcp(query, str+pivots[mid], str_end,
 285                                     lcp, maxlcp);
 286         if (newlcp == maxlcp) {
 287           lcp_left = lcp;
 288           lcp_right = minlcp;
 289           goto fullmatch;
 290         }
 291         if (query+newlcp == str_end ||
 292             query[newlcp] < str[pivots[mid]+newlcp]) {
 293           right = mid;
 294           minlcp = lcp;
 295           lcp_diff = newlcp - minlcp;
 296           if (lcp_diff == 0) goto in_the_middle;
 297           else goto closer_to_right;
 298         } else {
 299           left = mid+1;
 300           lcp_diff = minlcp - newlcp;
 301         }
 302       } else {
 303         left = mid+1;
 304       }
 305     }
 306   }
 307
 308  closer_to_right:
 309   while (true) {
 310     /*
 311     std::cout << "right minlcp=" << int(minlcp)
 312               << " lcp_diff=" << int(lcp_diff)
 313               << " mid=" << left + (right-left)/2
 314               << "\n";
 315     */
 316     if (left == right) return left;
 317     mid = left + (right-left)/2;
 318     if (str+pivots[mid] == query) return mid+1;
 319     lcptype lcp_mid = lcps[mid];
 320     if (lcp_mid < lcp_diff) {
 321       if (lcp_mid > 0) {
 322         lcp_diff -= lcp_mid;
 323         minlcp += lcp_mid;
 324       }
 325       left = mid+1;
 326     } else {
 327       if (lcp_mid == lcp_diff) {
 328         lcptype lcp = minlcp + lcp_diff;
 329         lcptype newlcp = suffix_lcp(query, str+pivots[mid], str_end,
 330                                     lcp, maxlcp);
 331         if (newlcp == maxlcp) {
 332           lcp_left = minlcp;
 333           lcp_right = lcp;
 334           goto fullmatch;
 335         }
 336         if (query+newlcp == str_end ||
 337             query[newlcp] < str[pivots[mid]+newlcp]) {
 338           right = mid;
 339           lcp_diff = newlcp - minlcp;
 340         } else {
 341           left = mid+1;
 342           minlcp = lcp;
 343           lcp_diff = minlcp - newlcp;
 344           if (lcp_diff == 0) goto in_the_middle;
 345           else goto closer_to_left;
 346         }
 347       }
 348       right = mid;
 349     }
 350   }
 351
 352  fullmatch:
 353   // find the range of fully matching pivots
 354   //std::cout << "fullmatch: left\n";
 355   {
 356     difftype right_tmp = mid;
 357     while (left != right_tmp) {
 358       difftype mid = left + (right_tmp-left)/2;
 359       //std::cout << left << mid << right_tmp << std::endl;
 360       difftype lcp_diff = lcps[mid];
 361       //std::cout << " " << lcp_left << " " << lcp_diff << std::endl;
 362       if (lcp_left + lcp_diff == maxlcp) {
 363         right_tmp = mid;
 364       } else {
 365         if (lcp_diff > 0) { lcp_left += lcp_diff; }
 366         left = mid+1;
 367       }
 368     }
 369   }
 370   //std::cout << "fullmatch: right\n";
 371   {
 372     difftype left_tmp = mid+1;
 373     while (left_tmp != right) {
 374       difftype mid = left_tmp + (right-left_tmp)/2;
 375       //std::cout << left_tmp << mid << right << std::endl;
 376       difftype lcp_diff = lcps[mid];
 377       //std::cout << " " << lcp_right << " " << lcp_diff << std::endl;
 378       if (lcp_right - lcp_diff == maxlcp) {
 379         left_tmp = mid+1;
 380       } else {
 381         if (lcp_diff < 0) { lcp_right -= lcp_diff; }
 382         right = mid;
 383       }
 384     }
 385   }
 386   //std::cout << "fullmatch: dcsearch\n";
 387   while (left != right) {
 388     difftype mid = left + (right-left)/2;
 389     //std::cout << left << mid << right << std::endl;
 390     if (dc_less_or_equal(pivots[mid], query-str)) {
 391       left = mid+1;
 392     } else {
 393       right = mid;
 394     }
 395   }
 396   return left;
 397
 398 }
 399
 400 #endif
 401
 402 };
 403
 404
 405 #endif // DISTRIBUTE_HPP