您的位置:首页 > 其它

codeforces 528D Fuzzy Search FFT

2016-12-13 21:00 441 查看
单调处理出第一个串往两边K个能覆盖的字符。

然后对四种字符每种字符做一遍FFT。

对于第一个串,如果这个位置不能匹配当前字符,那该位置为1。

对于第二个串,如果这个位置为当前字符,那该位置为1。

把第二个串反过来和第一个串跑FFT。

答案是在所有FFT中都为0的合法位置个数。

#include <bits/stdc++.h>
using namespace std;
#define N 210000
const double PI=acos(-1);
int n,m,K,len,ans;
char s
;
int a1
,b1
,cnt[5];
int val
[5],fin[N<<2];
int get(char c)
{return c=='A' ? 1:(c=='T' ? 2:(c=='G' ? 3:4));}
struct cp
{
double x,y;
cp(){}
cp(double x,double y):x(x),y(y){}
friend cp operator + (const cp &r1,const cp &r2)
{return cp(r1.x+r2.x,r1.y+r2.y);}
friend cp operator - (const cp &r1,const cp &r2)
{return cp(r1.x-r2.x,r1.y-r2.y);}
friend cp operator * (const cp &r1,const cp &r2)
{return cp(r1.x*r2.x-r1.y*r2.y,r1.x*r2.y+r2.x*r1.y);}
}a[N<<2],b[N<<2],c[N<<2];
void FFT(cp *a,int len,int type)
{
for(int i=0,t=0;i<len;i++)
{
if(i<t)swap(a[i],a[t]);
for(int j=len>>1;(t^=j)<j;j>>=1);
}
for(int i=2;i<=len;i<<=1)
{
cp wn(cos(2*PI*type/i),sin(2*PI*type/i));
for(int j=0;j<len;j+=i)
{
cp w(1,0),t;
for(int k=0;k<i>>1;k++,w=w*wn)
{
t=w*a[j+k+(i>>1)];
a[j+k+(i>>1)]=a[j+k]-t;
a[j+k]=a[j+k]+t;
}
}
}
if(type==-1)
for(int i=0;i<len;i++)c[i].x/=len;
}
int main()
{
//freopen("tt.in","r",stdin);
scanf("%d%d%d",&n,&m,&K);
scanf("%s",s+1);
for(int i=1;i<=n;i++)a1[i]=get(s[i]);
scanf("%s",s+1);
for(int i=1;i<=m;i++)b1[i]=get(s[i]);
for(int i=1,l=1,r=0;i<=n;i++)
{
while(l<i-K)cnt[a1[l]]--,l++;
while(r<i+K&&r<n)r++,cnt[a1[r]]++;
for(int j=1;j<=4;j++)
val[i][j]=cnt[j] ? 1:0;
}
for(len=1;len<n+m;len<<=1);
for(int now=1;now<=4;now++)
{
memset(a,0,sizeof(a));
memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)a[i-1].x=~val[i][now]&1;
for(int i=1;i<=m;i++)b[m-i].x=(b1[i]==now);
FFT(a,len,1);FFT(b,len,1);
for(int i=0;i<len;i++)c[i]=a[i]*b[i];
FFT(c,len,-1);
for(int i=0;i<len;i++)
fin[i]+=(int)(c[i].x+0.5);
}
for(int i=1;i+m-1<=n;i++)
if(!fin[m+i-2])ans++;
printf("%d\n",ans);
return 0;
}
内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: