题意:
DNA序列,在母串s中匹配模式串t,对于s中每个位置i,只要s[i-k]到s[i+k]中有c就认为匹配了c。求有多少个位置匹配了t。
分析:
这个字符串匹配的方式,什么kmp,各种自动机都不灵。
所以有一个邪门功夫,fft字符串匹配。(做过洛谷《残缺的字符串》一题的应该都不陌生,带通配符的匹配字符串可以用fft卷积来做)
首先,由于字符集大小只有4,所以我们可以对每个字符分别考虑。
根据题意,对于每个字符,我们预处理a[i]数组,代表i位置是否能匹配当前字符(左右k个能匹配也算)
之后b[i]数组,保存t[i]位置的字符是不是c。
然后将b数组倒过来,fft做一下卷积,把答案累计起来检查总和是不是m就可以判断是否匹配。
代码:
![](https://images.cnblogs.com/OutliningIndicators/ContractedBlock.gif)
1 #include<bits/stdc++.h> 2 #define db double 3 #define cp complex<db> 4 using namespace std; 5 const int N=530005; 6 const db pi=acos(-1); 7 int ans,ton[N],n,m,k,lm,L,r[N]; 8 char s[N],t[N];cp a[N],b[N]; 9 void init(){ 10 scanf("%d%d%d%s%s",&n,&m,&k,s,t); 11 } void fft(cp *f,int op){ 12 for(int i=0;i<lm;i++) 13 if(i>r[i]) swap(f[i],f[r[i]]); 14 for(int l=1;l<lm;l<<=1){ 15 cp wn(cos(pi/l),op*sin(pi/l)); 16 for(int i=0;i<lm;i+=(l<<1)){ 17 cp w(1,0); 18 for(int j=0;j<l;j++,w*=wn){ 19 cp T=w*f[i+j+l]; 20 f[i+j+l]=f[i+j]-T;f[i+j]+=T; 21 } 22 } 23 } if(op==-1) for(int i=0;i<lm;i++) f[i]/=lm; 24 } void calc(char c){ 25 memset(a,0,sizeof(a));memset(b,0,sizeof(b)); 26 int cnt=0; 27 for(int i=0;i<k;i++) cnt+=(s[i]==c); 28 for(int i=0;i<n;i++){ 29 if(i+k<n) cnt+=(s[i+k]==c); 30 if(i-k-1>=0) cnt-=(s[i-k-1]==c); 31 a[i]=cnt>0?1:0;cout<<a[i]<<" "; 32 } for(int i=0;i<m;i++) b[i]=(t[i]==c);puts(""); 33 fft(a,1);fft(b,1); 34 for(int i=0;i<lm;i++) a[i]*=b[i]; 35 fft(a,-1);for(int i=0;i<=n-m;i++) 36 ton[i]+=(int)(a[i+m-1].real()+0.5); 37 } void solve(){ 38 for(lm=1;lm<=n+m-2;lm<<=1,L++); 39 for(int i=0;i<lm;i++) 40 r[i]=(r[i>>1]>>1)|((i&1)<<(L-1)); 41 reverse(t,t+m);calc('A'); 42 calc('G');calc('C');calc('T'); 43 for(int i=0;i<=n-m;i++) if(ton[i]==m) ans++; 44 printf("%d ",ans); 45 } int main(){ 46 init();solve();return 0; 47 }