51nod1565 模糊搜索

来源:互联网 发布:网一网络加速器官网 编辑:程序博客网 时间:2024/06/06 16:49

对每一个字符单独考虑,O(n) 扫一边求出原串能匹配的位置,把模板串反序做fft就可以了。

#include<cstdio>#include<cmath>#include<cstring>#include<algorithm>using namespace std;const int maxn=1000000;const double pi=acos(-1);struct Complex{    double x,y;    Complex operator + (const Complex &c) const    {        return (Complex){x+c.x,y+c.y};    }    Complex operator - (const Complex &c) const    {        return (Complex){x-c.x,y-c.y};    }    Complex operator * (const Complex &c) const    {        return (Complex){x*c.x-y*c.y,x*c.y+y*c.x};    }    Complex operator / (const double &k) const    {        return (Complex){x/k,y/k};    }}f[maxn],g[maxn],w[maxn];char s[maxn],r[maxn],c[]={'A','G','C','T'};int ok[maxn],rev[maxn],n,m,k,l,t;void fft(Complex *a,int flag){    int x;    Complex t1,t2;    for (int i=0;i<l;i++)        if (rev[i]>i) swap(a[i],a[rev[i]]);    for (int i=0;i<t;i++)        for (int j=0;j<l;j+=1<<(i+1))        {            x=0;            for (int k=j;k<j+(1<<i);k++)            {                t1=a[k];                t2=a[k+(1<<i)]*w[x];                a[k]=t1+t2;                a[k+(1<<i)]=t1-t2;                x+=flag*(1<<(t-i-1));                if (x<0) x+=l;            }        }    if (flag==-1)        for (int i=0;i<l;i++)            a[i]=a[i]/l;}int main(){    //freopen("b.in","r",stdin);    int ans=0,now,num;    scanf("%d%d%d",&n,&m,&k);    scanf("%s",s);    scanf("%s",r);    for (int i=0;i<n;i++) ok[i]=1;    l=1,t=0;    while (l<n) l<<=1,t++;    for (int i=0;i<l;i++)        w[i]=(Complex){cos(2*pi*i/l),sin(2*pi*i/l)};    for (int i=0;i<l;i++)        for (int j=0;j<t;j++)            rev[i]|=((i>>j)&1)<<(t-j-1);    for (int x=0;x<4;x++)    {        for (int i=0;i<l;i++) f[i]=g[i]=(Complex){0,0};        now=0;        for (int i=0;i<n;i++)            if (s[i]==c[x])            {                f[i].x=1;                now=k;            }            else if (now)            {                f[i].x=1;                now--;            }        now=0;        for (int i=n-1;i>=0;i--)            if (s[i]==c[x])            {                f[i].x=1;                now=k;            }            else if (now)            {                f[i].x=1;                now--;            }        num=0;        for (int i=0;i<m;i++)            if (r[i]==c[x])                g[m-i-1].x=1,num++;        fft(f,1);        fft(g,1);        for (int i=0;i<l;i++) f[i]=f[i]*g[i];        fft(f,-1);        for (int i=0;i<=n-m;i++)            if ((int)(f[i+m-1].x+0.5)<num) ok[i]=0;    }    for (int i=0;i<=n-m;i++)        if (ok[i]) ans++;    printf("%d\n",ans);}
原创粉丝点击