51nod1565 模糊搜索

对每一个字符单独考虑, 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);
}
评论
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包
实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值