暑期学习之植入（l,d)模体发现PMSPc++（四）

再写了穷举之后，这两天思考着把PMSP写了，但是写的适配性不好，还只是局限于模体（9,2），以后还得改进适配性（黑脸）。

算法写的比较杂 先把伪代码 贴一下：

伪代码算法来自：霍红卫, 林帅, 于强,等. 基于MapReduce的模体发现算法[J]. 中国科技论文, 2012(7):487-494.

贴一下代码：

#include <iostream>
#include <string>
#include <stdlib.h>
#include <fstream>
#include<vector>
using namespace std;

int H_len(string a,string b,int l);//计算两个同长字符串的海明距离

int main()
{
    //读取txt数据文件
    ifstream ifile;
    ifile.open("data.txt");
    string line;
    string s[20];
    int panduan=0;
    while (getline(ifile, line)) {
            if(panduan%2==1){
                s[(panduan-1)/2]=line;
                //cout << s[(panduan-1)/2] << endl;
            }
        panduan = panduan + 1;
    }

    string Letter("ACGT");  //定义字母表
    int _len = 4;           //定义字母表长度
    int H_d = 2;            //d
    int LL = 9;             //l
    int n = 600;            //序列长度
    int t = 20;             //序列条数

    ofstream outf;
    outf.open("result.txt");

    //对第一行每一个l-mer做循环
    for(int a = 0;a < n - LL + 1;a++)
    {
        string t_moti = s[0].substr(a,LL);
        vector< vector<string> > l_mer_List;
        for(int b = 1;b < t;b++)
        {
            vector<string> l_mer0;
            //记录下跟tmoti海明距离小于2d的l-mer
            for(int c = 0;c < n - LL + 1;c++)
            {
                string cut = s[b].substr(c,LL);
                if(H_len(cut,t_moti,LL)<2*H_d || H_len(cut,t_moti,LL)==2*H_d)
                {
                    l_mer0.push_back(cut);
                }
            }
            l_mer_List.push_back(l_mer0);
        }
        /*int test_k = 0;
        for(int i = 0;i<l_mer_List.size();i++)
            for(int j = 0;j<l_mer_List[i].size();j++)
            {
                test_k = test_k + 1 ;
                cout<<l_mer_List[i][j]<<" "<<l_mer_List.size() <<endl;
            }
            */

        //对只有一个位置变化的l-mer做循环
        for(int i = 0;i < LL;i++)
        {
            string q_moti = t_moti;
            for(int j = 0;j<4;j++)
            {
                if(q_moti[i]!= Letter[j])
                {
                    int q_count = 0;
                    q_moti[i] = Letter[j];
                    for(int ck = 0;ck<l_mer_List.size();ck++)
                        for(int cr = 0;cr <l_mer_List[ck].size();cr++)
                        {
                            if(H_len(q_moti,l_mer_List[ck][cr],LL)<H_d || H_len(q_moti,l_mer_List[ck][cr],LL)==H_d)
                            {
                                q_count = q_count + 1;
                                break;
                            }
                        }
                    if(q_count==19)
                    {
                        cout<<q_moti<<"1"<<endl;
                        outf<<q_moti<<endl;
                        break;
                    }
                }
            }
        }
        //对有两个位置变化的l-mer做循环
        for(int i = 0;i < LL-1;i++)
        {
            for(int j = i+1;j<LL;j++)
            {
                  for(int ci = 0;ci<4;ci++)
                {
                     string q_moti = t_moti;
                    if(q_moti[i]!=Letter[ci])
                    {
                         for(int cj = 0;cj<4;cj++)
                         {
                             if(q_moti[i]!=Letter[cj])
                             {
                                int q_count = 0;
                                q_moti[i] = Letter[ci];
                                q_moti[j] = Letter[cj];
                                int sum = 0;
                                for(int ck = 0;ck<l_mer_List.size();ck++)
                                    for(int cr = 0;cr <l_mer_List[ck].size();cr++)
                                    {
                                        if(H_len(q_moti,l_mer_List[ck][cr],LL)<H_d || H_len(q_moti,l_mer_List[ck][cr],LL)==H_d)
                                        {
                                            q_count = q_count + 1;
                                            sum = H_len(q_moti,l_mer_List[ck][cr],LL) + sum ;
                                            break;
                                        }
                                    }
                                if(q_count==19)
                                {
                                    cout<<q_moti<<" "<<sum<<endl;
                                    outf<<q_moti<<endl;
                                    break;
                                }

                             }
                         }
                    }

                }
            }

        }


    }
    outf.close();
    return 0;
}
//计算两个同长字符串的海明距离
int H_len(string a,string b,int l)
{
    int temp_d = 0;
    for(int i = 0;i < l;i++)
    {
        if(a[i]!=b[i])
        {
            temp_d = temp_d + 1;
        }
    }
    return temp_d;
}

代码写的很冗余，也不规范，望各位大佬批评指正。