51nod 1565 模糊搜索

最新推荐文章于 2021-12-01 21:41:57 发布

原创最新推荐文章于 2021-12-01 21:41:57 发布 · 661 阅读

0 ·

CC 4.0 BY-SA版权

文章标签：

#C++ #FFT #思路

思路同时被 3 个专栏收录

216 篇文章

订阅专栏

其它OJ

61 篇文章

订阅专栏

FFT和它的小伙伴

10 篇文章

订阅专栏

本文介绍了一种使用快速傅立叶变换(FFT)进行模糊字符串匹配的方法，特别适用于基因序列匹配问题。通过将问题分解为针对每种字符类型的独立匹配过程，并利用FFT加速计算，该方法能在给定的模糊度范围内高效地找出子串在主串中的出现次数。

有两个基因串S和T，他们只包含AGCT四种字符。现在你要找出T在S中出现了几次。

有一个门限值k≥0。T在S的第i(1≤i≤|S|-|T|+1)个位置中出现的条件如下：把T的开头和S的第i个字符对齐，然后T中的每一个字符能够在S中找到一样的，且位置偏差不超过k的，那么就认为T在S的第i个位置中出现。也就是说对于所有的 j (1≤j≤|T|)，存在一个 p (1≤p≤|S|)，使得|(i+j-1)-p|≤k 和[p]=T[j]都成立。

例如，根据这样的定义"ACAT"出现在"AGCAATTCAT"的第2，3和6的位置。

如果k=0，那么这个就是经典的字符串匹配问题。

现在给定门限和两个基因串S，T，求出T在S中出现的次数。

Input

单组测试数据。
第一行有三个整数 |S|,|T|,k (1≤|T|≤|S|≤200000, 0≤k≤200000)，表示S的长度，T的长度，门限值。
第二行给出基因串S。
第三行给出基因串T。
两个串都只包含大写字母'A', 'T', 'G' 和'C'。

Output

输出一个整数，表示T在S中出现的次数。

Input示例

样例输入1
10 4 1
AGCAATTCAT
ACAT

Output示例

样例输出1
3

~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~

FFT+思路～

一共只有4种字符，所以我们把这四种分开匹配，再用&把所有答案综合起来。

因为是模糊匹配，所以在[x-k,x+k]范围内都是可以模糊地计算到x上的。所以我们把每个x的左右k位都划分成1，这样就能匹配到[x-k,x+k]的范围内了。

然后我们把b反向，用4次FFT，检验匹配是否可行。

这次写FFT我学了自己开struct表示complex的方法，如果用complex的话估计就过不了了～

数组要开大！

我的代码常数巨大，1s卡过……所以如果发现T掉了1个点，再交一次就能过辣。

#include<cstdio>
#include<cstring>
#include<iostream>
#include<cmath>
using namespace std;
#define pi acos(-1)

int n,m,n1,m1,r[800001],l,c[800001],id[128],A[4][800001],B[4][800001],ans[800001],totans,k;
char s[800001]; 

struct E{
	double a,b;
	E() {}
	E(double u,double v):a(u),b(v) {}
	E operator + (const E&u) const {return E(a+u.a,b+u.b);}
	E operator - (const E&u) const {return E(a-u.a,b-u.b);}
	E operator * (const E&u) const {return E(a*u.a-b*u.b,a*u.b+b*u.a);}
	E operator ! () {return E(a,-b);}
}a[800001],b[800001];

void fft(E *u,int v)
{
	for(int i=0;i<n;i++) if(r[i]>i) swap(u[r[i]],u[i]);
	for(int i=1;i<n;i<<=1)
	{
		E wn(cos(pi/i),sin(pi/i)*v);
		for(int j=0;j<n;j+=(i<<1))
		{
			E w(1,0);
			for(int k=0;k<i;k++,w=w*wn)
			{
				E x=u[j+k],y=u[i+j+k]*w;
				u[j+k]=x+y;u[i+j+k]=x-y;
			}
		}
	}
	if(v==-1) for(int i=0;i<n;i++) u[i].a/=n;
}

int main()
{
	id['A']=0;id['T']=1;id['C']=2;id['G']=3;
	scanf("%d%d%d",&n1,&m1,&k);
	scanf("%s",s);
	for(int i=0;i<n1;i++)
	{
		int kk=id[s[i]];
		A[kk][i]=1;
		for(int j=i+1;j<=min(n1-1,i+k);j++)
		  if(id[s[j]]==kk) break;
		  else A[kk][j]=1;
		for(int j=i-1;j>=max(0,i-k);j--)
		  if(A[kk][j]) break;
		  else A[kk][j]=1;
	}
	scanf("%s",s);
	for(int i=0;i<m1;i++) B[id[s[i]]][m1-i-1]=1;
	m=n1+m1;for(int i=m1-1;i<=m-2;i++) ans[i]=1;
	for(n=1;n<=m;n<<=1) l++;
	for(int i=0;i<n;i++) r[i]=(r[i>>1]>>1)|((i&1)<<(l-1));
	for(int k=0;k<4;k++)
	{
		memset(a,0,sizeof(a));
		memset(b,0,sizeof(b));
		for(int i=0;i<n1;i++) a[i]=E(A[k][i],0);
		for(int i=0;i<m1;i++) b[i]=E(B[k][i],0);
		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=0;i<=n;i++) c[i]=(int)(a[i].a+0.5);
		int now=0;
		for(int i=0;i<m1;i++) now+=B[k][i];
		for(int i=m1-1;i<=m-2;i++) ans[i]&=(c[i]==now);
	}
	for(int i=m1-1;i<=m-2;i++) totans+=ans[i];
	printf("%d\n",totans);
	return 0;
}