#include "pseud.h" void CountModuloCenter(Uint1Ptr simul_seq, int seqlen, int **MRC, short lmax) { int indel, i, j; short ok_flag; for(i=0;i<3;i++){ for(j=0;j=seqlen){ok_flag=0;break;} if(simul_seq[i+j]!=simul_seq[i+j+indel]){ok_flag=0;break;} } if(ok_flag){MRC[0][indel-1]++;} } } for(indel=-1;indel>=-lmax;indel--){ for(i=0;i=seqlen){ok_flag=0;break;} if(simul_seq[i+j]!=simul_seq[i+j-indel]){ok_flag=0;break;} } if(ok_flag){ for(j=0;j<-indel;j++){ if((i+j-2*indel)>=seqlen){ok_flag=0;break;} if(simul_seq[i+j-indel]!=simul_seq[i+j-2*indel]){ok_flag=0;break;} } if(ok_flag){MRC[2][-indel-1]++;} /*"CAGCAGCAG"*/ else{MRC[1][-indel-1]++;} /*"AACAGCAGTC"*/ } } } return; } void CountModuloCenterFilter(Uint1Ptr simul_seq, int seqlen, int **MRC, short lmax, int filter/*, FILE *fout*/) { int indel, i, j; short ok_flag; short *mrc_pos, *check; for(i=0;i<3;i++){ for(j=0;j=seqlen){ok_flag=0;break;} if(simul_seq[i+j]!=simul_seq[i+j+indel]){ok_flag=0;break;} } if(ok_flag){mrc_pos[i+indel] = 1;} } /*fprintf(fout,"No filter\nmrc_pos:"); for(i=0;i=2){check[i]=(short)(check[i]+1);check[i+indel]=(short)(check[i+indel]+1);} } for(i=0;i=1){mrc_pos[i]=1;}else{mrc_pos[i]=0;}} /*fprintf(fout,"filter 1\nmrc_pos:"); for(i=0;i=3){check[i]=(short)(check[i]+1);check[i+indel]=(short)(check[i+indel]+1);check[i+2*indel]=(short)(check[i+2*indel]+1);} } for(i=0;i=1){mrc_pos[i]=1;}else{mrc_pos[i]=0;}} /*fprintf(fout,"filter 2\nmrc_pos:"); for(i=0;i=2){check[i]=(short)(check[i]+1);check[i+indel]=(short)(check[i+indel]+1);} } for(i=0;i=1){mrc_pos[i]=1;}else{mrc_pos[i]=0;}} /*fprintf(fout,"filter 2\nmrc_pos:"); for(i=0;i=1){MRC[0][indel-1]++;} } /*deletion*/ for(i=0;i=1){check[i]=(short)(check[i]+1);check[i-indel]=(short)(check[i-indel]+1);} } for(i=0;i=2){MRC[2][indel-1]++;} if(check[i]==1){MRC[1][indel-1]++;} } } free(mrc_pos); free(check); return; } OptionsBlkPtr OptionsBlkFree(OptionsBlkPtr options) { int i, **MRC; if(options==NULL)return NULL; if(options->orig_seq){ options->orig_seq = MemFree(options->orig_seq); } if(options->target_seq){ options->target_seq = MemFree(options->target_seq); } if(options->pi){ options->pi = MemFree(options->pi); } if(options->PMat){ options->PMat = MemFree(options->PMat); } if(options->pID){ options->pID = MemFree(options->pID); } if(options->MRC){ MRC = options->MRC; for(i=0;i<3;i++){MRC[i]=MemFree(MRC[i]);} MRC = MemFree(MRC); } options = MemFree(options); return NULL; } OptionsBlkPtr GetOptions(char* infile, FILE *fout) { OptionsBlkPtr options; FILE *fin; char title[16]; Uint1Ptr orig_seq; Uint1Ptr target_seq; int x; int i, lmax2; double *pi, *PMat; double *pID, sum; int **MRC; if((fin = fopen(infile,"rt"))==NULL){ printf("Cannot open data file(press any key)\n");getch(); return NULL; } options = (OptionsBlkPtr)malloc(sizeof(OptionsBlk)); if(options==NULL){fclose(fin); return NULL;} fscanf(fin,"%[^=]=%ld, ",title,&(options->max_iter)); printf("%s=%ld\n",title,(options->max_iter)); fprintf(fout,"%s=%ld\n",title,(options->max_iter)); fscanf(fin,"%[^=]=%d, ",title,&(options->num_throw)); printf("%s=%d\n",title,(options->num_throw)); fprintf(fout,"%s=%d\n",title,(options->num_throw)); fscanf(fin,"%[^=]=%ld, ",title,&(options->max_idnum)); printf("%s=%ld\n",title,(options->max_idnum)); fprintf(fout,"%s=%ld\n",title,(options->max_idnum)); fscanf(fin,"%[^=]=%d, ",title,&x); options->lmax = (short)x; options->lmax0=options->lmax; printf("%s=%d\n",title,(options->lmax)); fprintf(fout,"%s=%d\n",title,(options->lmax)); fscanf(fin,"%[^=]=%d, ",title,&x); options->nmax = (short)x; printf("%s=%d\n",title,(options->nmax)); fprintf(fout,"%s=%d\n",title,(options->nmax)); fscanf(fin,"%[^=]=%d, ",title,&(options->select_num)); printf("%s=%d\n",title,(options->select_num)); fprintf(fout,"%s=%d\n",title,(options->select_num)); fscanf(fin,"%[^=]=%d, ",title,&(options->seqtype)); printf("%s=%d\n",title,(options->seqtype)); fprintf(fout,"%s=%d\n",title,(options->seqtype)); /*allocate orig_seq, target_seq*/ orig_seq = (Uint1Ptr)calloc(MaxSeqLen+1,sizeof(Uint1)); if(orig_seq==NULL){options = OptionsBlkFree(options); fclose(fin); return NULL;} options->orig_seq = orig_seq; target_seq = (Uint1Ptr)calloc(MaxSeqLen+1,sizeof(Uint1)); if(target_seq==NULL){options = OptionsBlkFree(options); fclose(fin); return NULL;} options->target_seq = target_seq; fscanf(fin,"%[^=]=%d, ",title,&(options->model)); printf("%s=%d\n",title,(options->model)); fprintf(fout,"%s=%d\n",title,(options->model)); fscanf(fin,"%[^=]=%d, ",title,&(options->num_t)); printf("%s=%d\n",title,(options->num_t)); fprintf(fout,"%s=%d\n",title,(options->num_t)); pi = (double *)malloc(NASize*sizeof(double)); for(i=0;ipi = pi; fscanf(fin,"%[^=]=%lf, ",title,&(options->kappa)); printf("%s=%1.2e\n",title,(options->kappa)); fprintf(fout,"%s=%1.2e\n",title,(options->kappa)); /*PMatij*/ PMat = (double *)calloc(NASize*NASize,sizeof(double)); if(PMat==NULL){options = OptionsBlkFree(options); fclose(fin); return NULL;} CalculatePMat(PMat,options,fout); options->PMat = PMat; /*pID*/ fscanf(fin,"%[^=]=%lf, ",title,&(options->lambda)); printf("%s=%1.2e\n",title,(options->lambda)); fprintf(fout,"%s=%1.2e\n",title,(options->lambda)); lmax2 = (int)(options->lmax)*2; pID = (double *)calloc(lmax2,sizeof(double)); if(pID==NULL){options = OptionsBlkFree(options); fclose(fin); return NULL;} /*for(i=0;ilmax);i++){ pID[i]=pow((double)(i+1), -1.7); sum += pID[i]; } for(i=0;i<(options->lmax);i++){ pID[i] /= sum; pID[(options->lmax)+i] = pID[i]*2.5; /*pseud gene pD = pI*2.5*/ } for(i=0;ipID = pID; fscanf(fin,"%lf, ",&(options->x)); printf("x=%1.2e\n",(options->x)); fprintf(fout,"z1=%1.2e\n",(options->x)); fscanf(fin,"%lf, ",&(options->y)); printf("y=%1.2e\n",(options->y)); fprintf(fout,"z2=%1.2e\n",(options->y)); fscanf(fin,"%lf, ",&(options->z)); printf("z=%1.2e\n",(options->z)); fprintf(fout,"z3=%1.2e\n",(options->z)); fscanf(fin,"%[^=]=%d, ",title, &(options->filter)); if((options->filter<0)||(options->filter>2)){ printf("Warning!Set filter within the limit(0-2)\n");getch(); options->filter=2; } printf("%s=%d\n",title, (options->filter)); fprintf(fout,"%s=%d\n",title, (options->filter)); MRC = (int **)calloc(3,sizeof(int *)); if(MRC==NULL){options = OptionsBlkFree(options); fclose(fin); return NULL;} for(i=0;i<3;i++){ MRC[i]=(int *)calloc(options->lmax,sizeof(int)); if(MRC[i]==NULL){options = OptionsBlkFree(options); fclose(fin); return NULL;} } options->MRC = MRC; fclose(fin); return options; }