【CF528D】Fuzzy Search（FFT）

【CF528D】Fuzzy Search（FFT）

题面

给定两个只含有(A,T,G,C)的(DNA)序列

定义一个字符(c)可以被匹配为：它对齐的字符，在距离(K)以内，存在一个字符(c)，问给定串(T)在(S)中出现了几次。

(|S|,|T|,K<=200000)

题解

字符集很小，可以分开进行(FFT)。

现在的匹配的定义为距离当前位置(K)以内的所有字符中是否含有这个字符，如果有设置为(1)，没有就是(0)，把字符分开做(FFT)然后相加，检查是否等于(|T|)即可。

#include<iostream>
#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<cmath>
#include<algorithm>
#include<set>
#include<map>
#include<vector>
#include<queue>
using namespace std;
#define ll long long
#define RG register
#define MAX 888888
const double Pi=acos(-1);
struct Complex{double a,b;}A[MAX],B[MAX],W[MAX];
Complex operator+(Complex a,Complex b){return (Complex){a.a+b.a,a.b+b.b};}
Complex operator-(Complex a,Complex b){return (Complex){a.a-b.a,a.b-b.b};}
Complex operator*(Complex a,Complex b){return (Complex){a.a*b.a-a.b*b.b,a.b*b.a+a.a*b.b};}
int r[MAX],N,n,m,l,K;
int ss[4][MAX],Ans[MAX];
char S[MAX],T[MAX],Box[4]={'A','T','G','C'};
void FFT(Complex *P,int opt)
{
	for(int i=1;i<N;++i)if(i<r[i])swap(P[i],P[r[i]]);
	for(int i=1;i<N;i<<=1)
		for(int p=i<<1,j=0;j<N;j+=p)
			for(int k=0;k<i;++k)
			{
				Complex w=(Complex){W[N/i*k].a,W[N/i*k].b*opt};
				Complex X=P[j+k],Y=w*P[i+j+k];
				P[j+k]=X+Y;P[i+j+k]=X-Y;
			}
	if(opt==-1)for(int i=0;i<N;++i)P[i].a/=N;
}
void Clear(){for(int i=0;i<N;++i)A[i].a=B[i].a=A[i].b=B[i].b=0;}
int main()
{
	scanf("%d%d%d",&n,&m,&K);
	for(N=1;N<=(n+m-2);N<<=1)++l;
	for(int i=0;i<N;++i)r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
	for(int i=1;i<N;i<<=1)
		for(int k=0;k<i;++k)W[N/i*k]=(Complex){cos(k*Pi/i),sin(k*Pi/i)};
	scanf("%s",S);scanf("%s",T);
	for(int i=1;i<=n;++i)
		for(int k=0;k<4;++k)
			if(S[i-1]==Box[k])ss[k][i]++;
	for(int i=1;i<=n;++i)
		for(int k=0;k<4;++k)ss[k][i]+=ss[k][i-1];
	for(int k=0;k<4;++k)
	{
		Clear();
		for(int i=0;i<n;++i)
			if(ss[k][min(n,i+K+1)]-ss[k][max(0,i-K)])
				A[i].a=1;
		for(int i=0;i<m;++i)
			if(T[m-i-1]==Box[k])B[i].a=1;
		FFT(A,1);FFT(B,1);
		for(int i=0;i<N;++i)A[i]=A[i]*B[i];
		FFT(A,-1);
		for(int i=m-1;i<n;++i)Ans[i-m+1]+=(int)(A[i].a+0.5);
	}
	int ans=0;
	for(int i=0;i<n;++i)if(Ans[i]==m)++ans;
	printf("%d
",ans);
	return 0;
}