interfac.c

clustalw1.83.DOS.ZIP,用于多序列比对的软件

           }
       else if (local_nseqs > 0)
           { 				/* success; found some seqs. */
		struct_penalties1 = NONE;
		if (sec_struct_mask1 != NULL) sec_struct_mask1=ckfree(sec_struct_mask1);
		if (gap_penalty_mask1 != NULL) gap_penalty_mask1=ckfree(gap_penalty_mask1);
		if (ss_name1 != NULL) ss_name1=ckfree(ss_name1);
                if (struct_penalties != NONE) /* feature table / mask in alignment */
                	{
					struct_penalties1 = struct_penalties;
					if (struct_penalties == SECST) {
						sec_struct_mask1 = (char *)ckalloc((max_aln_length) * sizeof (char));
						for (i=0;i<max_aln_length;i++)
							sec_struct_mask1[i] = sec_struct_mask[i];
					}
					gap_penalty_mask1 = (char *)ckalloc((max_aln_length) * sizeof (char));
					for (i=0;i<max_aln_length;i++)
						gap_penalty_mask1[i] = gap_penalty_mask[i];
        				ss_name1 = (char *)ckalloc( (MAXNAMES+1) * sizeof (char));

					strcpy(ss_name1,ss_name);
if (debug>0) {
for (i=0;i<seqlen_array[1];i++)
	fprintf(stdout,"%c",gap_penalty_mask1[i]);
fprintf(stdout,"\n");
}
					}
                nseqs = profile1_nseqs = local_nseqs;
				info("No. of seqs=%d",(pint)nseqs);
				profile1_empty=FALSE;
				profile2_empty=TRUE;
	   }
      }
    else
      {			        /* first seq to be read = profile1_nseqs + 1 */
       local_nseqs = readseqs(profile1_nseqs+(sint)1); 
       if(local_nseqs < 0)               /* file could not be opened */
           { 
		return local_nseqs;
           }
       else if(local_nseqs == 0)         /* no sequences */
           {
	       error("No sequences in file!  Bad format?");
		return local_nseqs;
           }
       else if(local_nseqs > 0)
           {
		struct_penalties2 = NONE;
		if (sec_struct_mask2 != NULL) sec_struct_mask2=ckfree(sec_struct_mask2);
		if (gap_penalty_mask2 != NULL) gap_penalty_mask2=ckfree(gap_penalty_mask2);
		if (ss_name2 != NULL) ss_name2=ckfree(ss_name2);
                if (struct_penalties != NONE) /* feature table / mask in alignment */
                	{
					struct_penalties2 = struct_penalties;
					if (struct_penalties == SECST) {
						sec_struct_mask2 = (char *)ckalloc((max_aln_length) * sizeof (char));
						for (i=0;i<max_aln_length;i++)
							sec_struct_mask2[i] = sec_struct_mask[i];
					}
					gap_penalty_mask2 = (char *)ckalloc((max_aln_length) * sizeof (char));
					for (i=0;i<max_aln_length;i++)
						gap_penalty_mask2[i] = gap_penalty_mask[i];
        				ss_name2 = (char *)ckalloc( (MAXNAMES+1) * sizeof (char));
					strcpy(ss_name2,ss_name);
if (debug>0) {
for (i=0;i<seqlen_array[profile1_nseqs+1];i++)
	fprintf(stdout,"%c",gap_penalty_mask2[i]);
fprintf(stdout,"\n");
}
					}
				info("No. of seqs in profile=%d",(pint)local_nseqs);
                nseqs = profile1_nseqs + local_nseqs;
                info("Total no. of seqs     =%d",(pint)nseqs);
				profile2_empty=FALSE;
				empty = FALSE;
	   }

      }
	if (sec_struct_mask != NULL) sec_struct_mask=ckfree(sec_struct_mask);
	if (gap_penalty_mask != NULL) gap_penalty_mask=ckfree(gap_penalty_mask);
	if (ss_name != NULL) ss_name=ckfree(ss_name);

	if(local_nseqs<=0) return local_nseqs;
	
	info("Sequences assumed to be %s",
		dnaflag?"DNA":"PROTEIN");
	if (usemenu) fprintf(stdout,"\n\n");
        for(i=profile2_empty?1:profile1_nseqs+1; i<=nseqs; i++) {
                info("Sequence %d: %-*s   %6.d %s",
                   (pint)i,max_names,names[i],(pint)seqlen_array[i],dnaflag?"bp":"aa");
        }	
	if(dnaflag) {
		gap_open   = dna_gap_open;
		gap_extend = dna_gap_extend;
	}
	else {
		gap_open   = prot_gap_open;
		gap_extend = prot_gap_extend;
	}

	return nseqs;
}



static void calc_gap_penalty_mask(int prf_length, char *mask, char *gap_mask)
{
	int i,j;
	char *struct_mask;

	struct_mask = (char *)ckalloc((prf_length+1) * sizeof(char));
/*
    calculate the gap penalty mask from the secondary structures
*/
	i=0;
	while (i<prf_length) {
		if (tolower(mask[i]) == 'a' || mask[i] == '$') {
			for (j = -helix_end_plus; j<0; j++) {
				if ((i+j>=0) && (tolower(struct_mask[i+j]) != 'a')
				             && (tolower(struct_mask[i+j]) != 'b'))
					struct_mask[i+j] = 'a';
			}
			for (j = 0; j<helix_end_minus; j++) {
				if (i+j>=prf_length || (tolower(mask[i+j]) != 'a'
				                    && mask[i+j] != '$')) break;
				struct_mask[i+j] = 'a';
			}
			i += j;
			while (tolower(mask[i]) == 'a'
				                    || mask[i] == '$') {
				if (i>=prf_length) break;
				if (mask[i] == '$') {
					struct_mask[i] = 'A';
					i++;
					break;
				}
				else struct_mask[i] = mask[i];
				i++;
			}
			for (j = 0; j<helix_end_minus; j++) {
				if ((i-j-1>=0) && (tolower(mask[i-j-1]) == 'a'
				                    || mask[i-j-1] == '$'))
					struct_mask[i-j-1] = 'a';
			}
			for (j = 0; j<helix_end_plus; j++) {
				if (i+j>=prf_length) break;
				struct_mask[i+j] = 'a';
			}
		}
	 	else if (tolower(mask[i]) == 'b' || mask[i] == '%') {
			for (j = -strand_end_plus; j<0; j++) {
				if ((i+j>=0) && (tolower(struct_mask[i+j]) != 'a')
				             && (tolower(struct_mask[i+j]) != 'b'))
					struct_mask[i+j] = 'b';
			}
			for (j = 0; j<strand_end_minus; j++) {
				if (i+j>=prf_length || (tolower(mask[i+j]) != 'b'
				                    && mask[i+j] != '%')) break;
				struct_mask[i+j] = 'b';
			}
			i += j;
			while (tolower(mask[i]) == 'b'
				                    || mask[i] == '%') {
				if (i>=prf_length) break;
				if (mask[i] == '%') {
					struct_mask[i] = 'B';
					i++;
					break;
				}
				else struct_mask[i] = mask[i];
				i++;
			}
			for (j = 0; j<strand_end_minus; j++) {
				if ((i-j-1>=0) && (tolower(mask[i-j-1]) == 'b'
				                    || mask[i-j-1] == '%'))
				struct_mask[i-j-1] = 'b';
			}
			for (j = 0; j<strand_end_plus; j++) {
				if (i+j>=prf_length) break;
 				struct_mask[i+j] = 'b';
			}
		}
	else i++;
	}

	for(i=0;i<prf_length;i++) {
		switch (struct_mask[i]) {
			case 'A':
				gap_mask[i] = helix_penalty+'0';
				break;
			case 'a':
				gap_mask[i] = helix_end_penalty+'0';
				break;
			case 'B':
				gap_mask[i] = strand_penalty+'0';
				break;
			case 'b':
				gap_mask[i] = strand_end_penalty+'0';
				break;
			default:
				gap_mask[i] = loop_penalty+'0';
				break;
		}
	}

	struct_mask=ckfree(struct_mask);
	
}

void print_sec_struct_mask(int prf_length, char *mask, char *struct_mask)
{
	int i,j;

/*
    calculate the gap penalty mask from the secondary structures
*/
	i=0;
	while (i<prf_length) {
		if (tolower(mask[i]) == 'a' || mask[i] == '$') {
			for (j = 0; j<helix_end_minus; j++) {
				if (i+j>=prf_length || (tolower(mask[i+j]) != 'a'
				                    && mask[i+j] != '$')) break;
				struct_mask[i+j] = 'a';
			}
			i += j;
			while (tolower(mask[i]) == 'a'
				                    || mask[i] == '$') {
				if (i>=prf_length) break;
				if (mask[i] == '$') {
					struct_mask[i] = 'A';
					i++;
					break;
				}
				else struct_mask[i] = mask[i];
				i++;
			}
			for (j = 0; j<helix_end_minus; j++) {
				if ((i-j-1>=0) && (tolower(mask[i-j-1]) == 'a'
				                    || mask[i-j-1] == '$'))
					struct_mask[i-j-1] = 'a';
			}
		}
	 	else if (tolower(mask[i]) == 'b' || mask[i] == '%') {
			for (j = 0; j<strand_end_minus; j++) {
				if (i+j>=prf_length || (tolower(mask[i+j]) != 'b'
				                    && mask[i+j] != '%')) break;
				struct_mask[i+j] = 'b';
			}
			i += j;
			while (tolower(mask[i]) == 'b'
				                    || mask[i] == '%') {
				if (i>=prf_length) break;
				if (mask[i] == '%') {
					struct_mask[i] = 'B';
					i++;
					break;
				}
				else struct_mask[i] = mask[i];
				i++;
			}
			for (j = 0; j<strand_end_minus; j++) {
				if ((i-j-1>=0) && (tolower(mask[i-j-1]) == 'b'
				                    || mask[i-j-1] == '%'))
				struct_mask[i-j-1] = 'b';
			}
		}
	else i++;
	}
}



FILE *  open_output_file(char *prompt,      char *path, 
				char *file_name,   char *file_extension)
 
{	static char temp[FILENAMELEN+1];
	static char local_prompt[MAXLINE];
	FILE * file_handle;

/*	if (*file_name == EOS) {
*/		strcpy(file_name,path);
		strcat(file_name,file_extension);
/*	}
*/
	if(strcmp(file_name,seqname)==0) {
		warning("Output file name is the same as input file.");
		if (usemenu) {
			strcpy(local_prompt,"\n\nEnter new name to avoid overwriting ");
			strcat(local_prompt," [%s]: ");          
			fprintf(stdout,local_prompt,file_name);
			gets(temp);
			if(*temp != EOS) strcpy(file_name,temp);
		}
	}
	else if (usemenu) {
		strcpy(local_prompt,prompt);
		strcat(local_prompt," [%s]: ");          
		fprintf(stdout,local_prompt,file_name);
		gets(temp);
		if(*temp != EOS) strcpy(file_name,temp);
	}

#ifdef VMS
	if((file_handle=fopen(file_name,"w","rat=cr","rfm=var"))==NULL) {
#else
	if((file_handle=fopen(file_name,"w"))==NULL) {
#endif
		error("Cannot open output file [%s]",file_name);
		return NULL;
	}
	return file_handle;
}



FILE *  open_explicit_file(char *file_name)
{ 
	FILE * file_handle;

	if (*file_name == EOS) {
		error("Bad output file [%s]",file_name);
		return NULL;
	}
#ifdef VMS
	if((file_handle=fopen(file_name,"w","rat=cr","rfm=var"))==NULL) {
#else
	if((file_handle=fopen(file_name,"w"))==NULL) {
#endif
		error("Cannot open output file [%s]",file_name);
		return NULL;
	}
	return file_handle;
}


/* Ramu void */

void align(char *phylip_name)
{ 
	char path[FILENAMELEN+1];
	FILE *tree;
	sint count;
	
	if(empty && usemenu) {
		error("No sequences in memory. Load sequences first.");
		return;
	}

	   struct_penalties1 = struct_penalties2 = NONE;
	   if (sec_struct_mask1 != NULL) sec_struct_mask1=ckfree(sec_struct_mask1);
	   if (sec_struct_mask2 != NULL) sec_struct_mask2=ckfree(sec_struct_mask2);
	   if (gap_penalty_mask1 != NULL) gap_penalty_mask1=ckfree(gap_penalty_mask1);
	   if (gap_penalty_mask2 != NULL) gap_penalty_mask2=ckfree(gap_penalty_mask2);
	   if (ss_name1 != NULL) ss_name1=ckfree(ss_name1);
	   if (ss_name2 != NULL) ss_name2=ckfree(ss_name2);


        get_path(seqname,path);
/* DES DEBUG 
	fprintf(stdout,"\n\n Seqname = %s  \n Path = %s \n\n",seqname,path);
*/
	if(usemenu || !interactive) {
        	if(!open_alignment_output(path)) return;
	}

	if (nseqs >= 2) {

        	get_path(seqname,path);
        	if (phylip_name[0]!=EOS) {
                	if((tree = open_explicit_file(
                	phylip_name))==NULL) return;
        	}
        	else {
                 	if((tree = open_output_file(
                	"\nEnter name for new GUIDE TREE           file  ",path,
                	phylip_name,"dnd")) == NULL) return;
        	}
	}

	if (save_parameters) create_parameter_output();

	if(reset_alignments_new || reset_alignments_all) reset_align();

        info("Start of Pairwise alignments");
        info("Aligning...");
        if(dnaflag) {
                gap_open   = dna_gap_open;
                gap_extend = dna_gap_extend;
                pw_go_penalty  = dna_pw_go_penalty;
                pw_ge_penalty  = dna_pw_ge_penalty;
                ktup       = dna_ktup;
                window     = dna_window;
                signif     = dna_signif;
                wind_gap   = dna_wind_gap;

        }
        else {
                gap_open   = prot_gap_open;
                gap_extend = prot_gap_extend;
                pw_go_penalty  = prot_pw_go_penalty;
                pw_ge_penalty  = prot_pw_ge_penalty;
                ktup       = prot_ktup;
                window     = prot_window;
                signif     = prot_signif;
                wind_gap   = prot_wind_gap;

        }

        if (quick_pairalign)
           show_pair((sint)0,nseqs,(sint)0,nseqs);
        else
           pairalign((sint)0,nseqs,(sint)0,nseqs);

	if (nseqs >= 2) {

		guide_tree(tree,1,nseqs);
		info("Guide tree        file created:   [%s]",
                phylip_name);
	}

	
	count = malign((sint)0,phylip_name);
	
	if (count <= 0) return;

	if (usemenu) fprintf(stdout,"\n\n\n");
	
	create_alignment_output(1,nseqs);
        if (showaln && usemenu) show_aln();
	phylip_name[0]=EOS;
	return ;
}





void new_sequence_align(char *phylip_name)
{ 
	char path[FILENAMELEN+1];
	char tree_name[FILENAMELEN+1],temp[MAXLINE+1];
	Boolean use_tree;
	FILE *tree;
	sint i,j,count;
	float dscore;
	Boolean save_ss2;
	
	if(profile1_empty && usemenu) {
		error("No profile in memory. Input 1st profile first.");
		return;
	}

	if(profile2_empty && usemenu) {
		error("No sequences in memory. Input sequences first.");
		return;
	}

        get_path(profile2_name,path);

        if(usemenu || !interactive) {
        	if(!open_alignment_output(path)) return;
	}

	new_seq = profile1_nseqs+1;

/* check for secondary structure information for list of sequences */

	save_ss2 = use_ss2;
	if (struct_penalties2 != NONE &&