| 308 | } |
| 309 | |
| 310 | static Support |
| 311 | testSequence(const Sequence& sequence) |
| 312 | { |
| 313 | static unsigned char BASES[] = { 'A', 'C', 'T', 'G' }; |
| 314 | int found = 0; |
| 315 | int tests = 0; |
| 316 | unsigned r = g_vanillaBloom->get_k(); |
| 317 | if (sequence.size() >= r) { |
| 318 | tests = sequence.size() - r + 1; |
| 319 | int offset = 0; |
| 320 | if (opt::errorCorrection) { |
| 321 | btllib::NtHash nthash(sequence, g_vanillaBloom->get_hash_num(), r); |
| 322 | for (const auto& hitSeeds : g_spacedSeedsBloom->contains(sequence)) { |
| 323 | nthash.roll(); |
| 324 | if (hitSeeds.size() > 0) { |
| 325 | nthash.sub({}, {}); |
| 326 | if (g_vanillaBloom->contains(nthash.hashes())) { |
| 327 | found++; |
| 328 | } else { |
| 329 | bool success = false; |
| 330 | for (const auto& hitSeed : hitSeeds) { |
| 331 | const auto seed = g_spacedSeedsBloom->get_parsed_seeds()[hitSeed]; |
| 332 | for (auto seedIt = |
| 333 | (seed.begin() + |
| 334 | std::round( |
| 335 | seed.size() * (1.00 - SPACED_SEEDS_SNP_FRACTION))); |
| 336 | seedIt != seed.end(); |
| 337 | ++seedIt) { |
| 338 | const auto pos = *seedIt; |
| 339 | for (auto base : BASES) { |
| 340 | if (base == (unsigned char)(sequence[offset + pos])) { |
| 341 | continue; |
| 342 | } |
| 343 | nthash.sub({ pos }, { base }); |
| 344 | if (g_vanillaBloom->contains(nthash.hashes())) { |
| 345 | success = true; |
| 346 | found++; |
| 347 | break; |
| 348 | } |
| 349 | } |
| 350 | if (success) { |
| 351 | break; |
| 352 | } |
| 353 | } |
| 354 | if (success) { |
| 355 | break; |
| 356 | } |
| 357 | } |
| 358 | } |
| 359 | } |
| 360 | offset++; |
| 361 | } |
| 362 | } else { |
| 363 | found = g_vanillaBloom->contains(sequence); |
| 364 | } |
| 365 | } |
| 366 | return Support(found, tests); |
| 367 | } |