B - Fuzzy Search （FFT）

  1 #include<iostream>
  2 #include<cstring>
  3 #include<string>
  4 #include<cmath>
  5 #include<algorithm>
  6 #include<stdio.h>
  7 using namespace std;
  8 # define ll long long
  9 const double PI = acos(-1.0);
 10 const int maxn = 8e5+100;
 11 struct complex
 12 {
 13     double r,i;
 14     complex(double _r = 0,double _i = 0)
 15     {
 16         r = _r;
 17         i = _i;
 18     }
 19     complex operator +(const complex &b)
 20     {
 21         return complex(r+b.r,i+b.i);
 22     }
 23     complex operator -(const complex &b)
 24     {
 25         return complex(r-b.r,i-b.i);
 26     }
 27     complex operator *(const complex &b)
 28     {
 29         return complex(r*b.r-i*b.i,r*b.i+i*b.r);
 30     }
 31 };
 32 void change(complex y[],int len)
 33 {
 34     int i,j,k;
 35     for(i = 1, j = len/2; i < len-1; i++)
 36     {
 37         if(i < j)
 38             swap(y[i],y[j]);
 39         k = len/2;
 40         while( j >= k)
 41         {
 42             j -= k;
 43             k /= 2;
 44         }
 45         if(j < k)
 46             j += k;
 47     }
 48 }
 49 void fft(complex y[],int len,int on)
 50 {
 51     change(y,len);
 52     for(int h = 2; h <= len; h <<= 1)
 53     {
 54         complex wn(cos(-on*2*PI/h),sin(-on*2*PI/h));
 55         for(int j = 0; j < len; j += h)
 56         {
 57             complex w(1,0);
 58             for(int k = j; k < j+h/2; k++)
 59             {
 60                 complex u = y[k];
 61                 complex t = w*y[k+h/2];
 62                 y[k] = u+t;
 63                 y[k+h/2] = u-t;
 64                 w = w*wn;
 65             }
 66         }
 67     }
 68     if(on == -1)
 69         for(int i = 0; i < len; i++)
 70             y[i].r /= len;
 71 }
 72 char str1[maxn],str2[maxn];
 73 complex x1[maxn],x2[maxn];
 74 int len1,len2,len=1;
 75 int vis[maxn][5],id[maxn],cnt[10];
 76 ll ans[maxn];
 77 int main()
 78 {
 79     int n,m,k;
 80     scanf("%d %d %d",&n,&m,&k);
 81     scanf("%s",str1+1);
 82     scanf("%s",str2+1);
 83     len1=strlen(str1+1);
 84     len2=strlen(str2+1);
 85     while(len<len1*2||len<len2*2)
 86         len<<=1;
 87    //     cout<<len<<endl;
 88     int l=0,r=0;
 89     id[‘A‘]=1,id[‘G‘]=2,id[‘T‘]=3,id[‘C‘]=4;
 90     for(int i=1; i<=n; i++)
 91     {
 92         while(l<n&&l<i-k)//判断这个区间内的字符。
 93             cnt[id[(int)str1[l++]]]--;
 94         while(r<n&&r<i+k)
 95             cnt[id[(int)str1[++r]]]++;
 96         for(int j=1; j<=4; j++)
 97             if(cnt[j])
 98                 vis[i][j]=1;
 99     }
100     for(int i=1; i<=4; i++)
101     {
102         for(int j=0; j<len; j++)
103         {
104             x1[j]=complex(0,0);
105             x2[j]=complex(0,0);
106         }
107         for(int j=1; j<=n; j++)
108         {
109             if(vis[j][i])
110                 x1[j-1]=complex(1,0);
111         }
112         for(int j=1; j<=m; j++)
113         {
114             if(id[(int)str2[j]]==i)
115                 x2[m-j]=complex(1,0);
116         }
117         fft(x1,len,1);
118         fft(x2,len,1);
119         for(int j=0; j<len; j++)
120         {
121             x1[j]=x1[j]*x2[j];
122         }
123         fft(x1,len,-1);
124         for(int j=0; j<len; j++)
125         {
126             ans[j]+=(ll)(x1[j].r+0.5);
127         }
128     }
129     int num=0;
130     for(int i=0; i<len; i++)
131     {
132     //    cout<<i<<" "<<ans[i]<<endl;
133         if(ans[i]==m)
134             num++;
135     }
136     printf("%d\n",num);
137     return 0;
138 }