Advertisement
Guest User

Untitled

a guest
Apr 20th, 2019
85
0
Never
Not a member of Pastebin yet? Sign Up, it unlocks many cool features!
text 2.97 KB | None | 0 0
  1. #include <iostream>
  2. #include <vector>
  3. #include <algorithm>
  4. #include <complex>
  5. using namespace std;
  6.  
  7. const long double M_PI = 3.141592653589793238462643383279;
  8.  
  9. struct DNA {
  10. DNA(int M) : A(M, 0), C(M, 0), G(M, 0), T(M, 0) {}
  11.  
  12. vector<complex<long double>> A, C, G, T;
  13. };
  14.  
  15. void dfft(vector<complex<long double>>& vec, bool inv) {
  16. int sz = vec.size(), bit;
  17. for (int i = 1, j = 0; i < sz; i++) {
  18. bit = sz >> 1;
  19. for (; j >= bit; bit >>= 1)
  20. j -= bit;
  21. j += bit;
  22. if (i < j)
  23. swap(vec[i], vec[j]);
  24. }
  25. for (int len = 2; len <= sz; len <<= 1) {
  26. long double ang = 2 * M_PI / len * (inv ? -1 : 1);
  27. complex<long double> wlen(cos(ang), sin(ang));
  28. for (int i = 0; i < sz; i += len) {
  29. complex<long double> w(1);
  30. for (int j = 0; j < len >> 1; j++) {
  31. complex<long double> u = vec[i + j], v = w * vec[i + j + (len >> 1)];
  32. vec[i + j] = u + v;
  33. vec[i + j + (len >> 1)] = u - v;
  34. w *= wlen;
  35. }
  36. }
  37. }
  38. if (inv)
  39. for (int i = 0; i < sz; i++)
  40. vec[i] /= sz;
  41. }
  42.  
  43. void parse(const string & s, DNA & d) {
  44. int t = 0;
  45.  
  46. for (auto it = s.begin(); it != s.end(); it++, t++)
  47. if (*it == 'A')
  48. d.A[t] = 1;
  49. else if (*it == 'G')
  50. d.G[t] = 1;
  51. else if (*it == 'C')
  52. d.C[t] = 1;
  53. else if (*it == 'T')
  54. d.T[t] = 1;
  55. }
  56.  
  57.  
  58. int main() {
  59. ios_base::sync_with_stdio(0);
  60. cin.tie(0);
  61. cout.tie(0);
  62.  
  63. int M;
  64. string s1, s2;
  65. cin >> M >> s1 >> s2;
  66. reverse(s1.begin(), s1.end());
  67. DNA first(M << 1), second(M << 1);
  68. parse(s1, first), parse(s2, second);
  69.  
  70. vector<complex<long double>> a_mul(M << 1, 0), c_mul(M << 1, 0), g_mul(M << 1, 0), t_mul(M << 1, 0);
  71. vector<unsigned long long int> a_real(M << 1, 0), c_real(M << 1, 0), g_real(M << 1, 0), t_real(M << 1, 0);
  72.  
  73. vector<unsigned long long> result(M, 0);
  74.  
  75. for (int i = 0; i < M; i++) {
  76. second.A[M + i] = second.A[i];
  77. second.C[M + i] = second.C[i];
  78. second.G[M + i] = second.G[i];
  79. second.T[M + i] = second.T[i];
  80. }
  81.  
  82. dfft(first.A, false), dfft(second.A, false);
  83. dfft(first.C, false), dfft(second.C, false);
  84. dfft(first.G, false), dfft(second.G, false);
  85. dfft(first.T, false), dfft(second.T, false);
  86.  
  87. for (int i = 0; i < M << 1; i++) {
  88. a_mul[i] = first.A[i] * second.A[i];
  89. c_mul[i] = first.C[i] * second.C[i];
  90. g_mul[i] = first.G[i] * second.G[i];
  91. t_mul[i] = first.T[i] * second.T[i];
  92. }
  93.  
  94. dfft(a_mul, true); dfft(c_mul, true);
  95. dfft(g_mul, true); dfft(t_mul, true);
  96.  
  97. for (int i = 0; i < M << 1; i++) {
  98. a_real[i] = (unsigned long long int) (a_mul[i].real() + 0.5);
  99. c_real[i] = (unsigned long long int) (c_mul[i].real() + 0.5);
  100. g_real[i] = (unsigned long long int) (g_mul[i].real() + 0.5);
  101. t_real[i] = (unsigned long long int) (t_mul[i].real() + 0.5);
  102. }
  103.  
  104. for (int i = M; i < M << 1; i++) {
  105. unsigned long long int t = a_real[i] + c_real[i] + g_real[i] + t_real[i];
  106. result[i - M] = t;
  107. }
  108.  
  109. cout << *max_element(result.begin(), result.end()) << ' ';
  110. cout << ((result.end() - max_element(result.begin(), result.end()) - 1) % M);
  111. return 0;
  112. }
Advertisement
Add Comment
Please, Sign In to add comment
Advertisement