/* Last changed Time-stamp: <2007-12-05 13:55:42 ronny> */ /* Ineractive Access to folding Routines c Ivo L Hofacker Vienna RNA package */ #include #include #include #include #include #include #include #include "snofold.h" #include "fold.h" #include "energy_par.h" #include "snoop.h" #include "part_func.h" #include "fold_vars.h" #include "PS_dot.h" #include "utils.h" #include "aln_util.h" #include "RNAsnoop_cmdl.h" static void aliprint_struc(snoopT *dup, const char **s1, const char **s2, char**,char**,int count, int); static void print_struc(snoopT *dup, const char *s1, const char *s2, char*, char*, int,int); /* static void print_struc_L(snoopT const *dup, const char *s1, const char *s2, char*, char*); */ static void redraw_output(char *fname, char *output, int plfold_up_flag, char *suffix, char *access, char *sname, char *tname); static int ** read_plfold_i(char *fname, const int beg, const int end); /* read plfold profiles */ static int ** read_rnaup(char *fname, const int beg, const int end); static int get_max_u(const char *s, char delim); extern int cut_point; /*@unused@*/ static char UNUSED rcsid[] = "$Id: RNAfold.c,v 1.23 2007/12/05 13:04:08 ivo Exp $"; #define PRIVATE static #define MAX_NUM_NAMES 500 /*--------------------------------------------------------------------------*/ int main(int argc, char *argv[]) { struct RNAsnoop_args_info args_info; char *string_s, *line_s, *name_s, *temp_s; /*string for snoRNA*/ char *string_t, *line_t, *temp_t; /*string for target RNA*/ char *structure=NULL, *cstruc=NULL; char *sname=NULL,*tname=NULL/*name of the sno RNa file, mRNA file respectively*/; char *name=NULL; char *access;char *result_file; char *output_directory; int fullStemEnergy; output_directory=NULL; result_file=NULL; access=NULL; char *suffix; suffix=NULL; FILE *sno, *mrna; int fast;fast=0; /* long int elapTicks; */ /* clock_t Begin, End; */ int plfold_up_flag=0; int nice,i,r,l,length_s,/* , length_t, */ penalty, /*extension penalty*/ threshloop, /*energy threshold on loop*/ threshLE, /*energy threshold on the S2 part*/ threshRE, /*energy threshold on the S1 part*/ threshDE, /*Total duplex energy threshold*/ threshTE, /*Total duplex + Loop energy*/ threshSE, /*Sum of all energies*/ threshD, /*lower stem energy*/ half_stem, /*minimal length of stem*/ max_half_stem, /*maximal length of stem*/ max_s2 /*maximale position wo stem anfangen darf in 3->5*/, min_s2 /*minimale position wo stem anfangen darf in 3->5*/, max_s1 /*minimal position wo stem enden darf in 3->5 */, min_s1 /*maximale position wo stem enden darf in 3->5*/, distance /*distance between two subopts*/, min_d1 /*minimal distance between 5' sno and first duplex*/, min_d2 /*minimal distance between 3' sno and second duplex*/, max_asymm /*maximal asymmetry in the stem interior loop*/, alignment_length/*maximal target RNA alignment length*/, alignment /*flag to use ali version*/, delta, redraw /*if used (option I) allow to redraw command line output into ps files */; int noconv=0; string_s=NULL; string_t=NULL; plfold_up_flag=0; alignment=0; redraw=0; nice=0; fast=0; snoop_subopt_sorted=0; /* ############################################# # check the command line parameters ############################################# */ if(RNAsnoop_cmdline_parser(argc,argv, &args_info)!=0) exit(1); /* Redraw results from RNAsnoop */ if(args_info.constraint_given) fold_constrained=1; /* structure constraint */ if(args_info.produce_ps_given) redraw=1; /*output directory for RNAsnoop results */ output_directory=strdup(args_info.output_directory_arg); /* Draw directly nice pictures from RNAsnoop*/ if(args_info.direct_redraw_given) nice=1; /*Accessibility files are coming from RNAup*/ if(args_info.from_RNAup_given) {plfold_up_flag=2;access=strdup(args_info.from_RNAup_arg);} /*Suffix for the name of the accessibility files*/ suffix=strdup(args_info.suffix_arg); /*Accessibility files are coming from RNAup*/ if(args_info.from_RNAplfold_given) {plfold_up_flag=1;access=strdup(args_info.from_RNAplfold_arg);} /*Check if we are in alignment mode*/ if(args_info.alignment_mode_given) alignment=1; /*Name of the stems file*/ if(args_info.query_given) sname=strdup(args_info.query_arg); /*Name of the target file*/ if(args_info.target_given) tname=strdup(args_info.target_arg); /*Delta for suboptimals*/ delta=(int) (100*(args_info.energy_threshold_arg+0.1)); /*fast folding and target search*/ fast=args_info.fast_folding_arg; /*Penalty for duplex extension*/ penalty=args_info.extension_cost_arg; /*threshold on loop energy*/ threshloop=args_info.minimal_loop_energy_arg; /*threshold on right duplex*/ threshRE=args_info.minimal_right_duplex_arg; /*threshold on left duplex*/ threshLE=args_info.minimal_left_duplex_arg; /*threshold on minimal duplex*/ threshDE=args_info.minimal_duplex_arg; /*threshold on duplex distance*/ distance=args_info.duplex_distance_arg; /*threshold on minimal stem length*/ half_stem=args_info.minimal_stem_length_arg; /*threshold on maximal stem length*/ max_half_stem=args_info.maximal_stem_length_arg; /*threshold on minimal duplex box length*/ min_s2=args_info.minimal_duplex_box_length_arg; /*threshold on maximal duplex box length*/ max_s2=args_info.maximal_duplex_box_length_arg; /*threshold on minimal snoRNA stem loop length*/ min_s1=args_info.minimal_snoRNA_stem_loop_length_arg; /*thrshold on maximal snoRNA stem loop length*/ max_s1=args_info.maximal_snoRNA_stem_loop_length_arg; /*threshold on minimal snoRNA duplex length*/ min_d1=args_info.minimal_snoRNA_duplex_length_arg; /*threshold on minimal snoRNA duplex length*/ min_d2=args_info.minimal_snoRNA_duplex_length_arg; /*threshold on minimal duplex stem energy*/ threshTE=args_info.minimal_duplex_stem_energy_arg; /*threshold on minimal total energy*/ threshSE=args_info.minimal_total_energy_arg; /*threshold on maximal stem asymmetry*/ max_asymm=args_info.maximal_stem_asymmetry_arg; /*threshold on minimal lower stem energy*/ threshD=args_info.minimal_lower_stem_energy_arg; /*threshold on minimal lower stem energy*/ alignment_length=args_info.alignmentLength_arg; threshloop=MIN2(threshloop,0); /* if(plfold_up_flag && !fast){ */ /* printf("Sorry, no accessibility implementation with the non fast implementation\n"); */ /* printf("If you want to include accessibility information please run RNAsnoop with the -f option\n"); */ /* printf("If you want to run RNAsnoop with the slow algorithm, please remove run RNAsnoop without -P\n"); */ /* exit(0); */ /* } */ if(plfold_up_flag==2 && suffix==NULL){ printf("RNAsnoop needs a suffix (-S u1-to-30) for the RNAup accessibility file in order to localize them\n"); exit(0); } if(redraw){ if(!alignment){ redraw_output(result_file, output_directory, plfold_up_flag, suffix, access,NULL,NULL); } else{ redraw_output(result_file, output_directory, plfold_up_flag, suffix, access,sname,tname); } /* readfile */ /* parse lines */ /* use current ploting function to output the results */ exit(0); } /* istty = isatty(fileno(stdout))&&isatty(fileno(stdin)); */ min_s2+=5; max_s2+=5; min_s1+=5; max_s1+=5; min_d1+=5; min_d2+=5; if(!alignment) { if(tname==NULL || sname==NULL) RNAsnoop_cmdline_parser_print_help(); sno=fopen(sname, "r"); if(sno==NULL){printf("%s: Wrong snoRNA file name\n", sname);return 0;} mrna=fopen(tname, "r"); if(mrna==NULL){printf("%s: Wrong target file name\n", tname);return 0;} do { /* main loop: continue until end of file */ if ((line_s = get_line(sno))==NULL) { free(line_s); break; } /* skip comment lines and get filenames */ while ((*line_s=='*')||(*line_s=='\0')||(*line_s=='>')) { if (*line_s=='>') name_s = (char *) space(strlen(line_s)+1); (void) sscanf(line_s,"%s",name_s); free(line_s); if ((line_s = get_line(sno))==NULL) { free(line_s); break; } } if(name_s == NULL) {printf("Your snoRNA sequence: \n%s\nhas no header. Please update your fasta file\n", line_s); exit(0);} /* if ((line ==NULL) || (strcmp(line, "@") == 0)) break; */ temp_s = (char *) space(strlen(line_s)+1); (void) sscanf(line_s,"%s",temp_s); free(line_s); length_s = (int) strlen(temp_s); for (l = 0; l < length_s; l++) { temp_s[l] = toupper(temp_s[l]); if (!noconv && temp_s[l] == 'T') temp_s[l] = 'U'; } string_s = (char*) space(length_s + 11); strcpy(string_s, "NNNNN"); strcat(string_s+5, temp_s); strcat(string_s+5+length_s, "NNNNN\0"); free(temp_s); /* We declare the structure variable here as it will also contains the final stem structure */ structure=(char *)space((unsigned) length_s+11); if(fold_constrained){ cstruc = get_line(sno); if(cstruc!=NULL){ int dn3=strlen(cstruc)-(length_s-10); strcpy(structure,"....."); strcat(structure, cstruc); if(dn3>=0){ strcat(structure,".....\0"); } else{ while(dn3++){ strcat(structure,"."); } strcat(structure,"\0"); } /* Now we fold with constraints the */ } } fullStemEnergy = snofold(string_s, structure, max_asymm, threshloop, min_s2, max_s2, half_stem, max_half_stem); do{ /* main loop for target continue until end of file */ snoopT mfe; if ((line_t = get_line(mrna))==NULL) { /* free(line_t); */ break; } char *name_t=NULL; /* skip comment lines and get filenames */ while ((*line_t=='*')||(*line_t=='\0')||(*line_t=='>')) { if (*line_t=='>'){ printf("%s\n", name_s); name_t = (char*) space(strlen(line_t)+1); (void) sscanf(line_t,"%s",name_t); printf("%s\n", name_t); free(line_t); /* free(string_t); */ } if ((line_t = get_line(mrna))==NULL) { free(line_t); break; } } if(name_t == NULL) {printf("Your target sequence: \n%s\nhas no header. Please update your fasta file\n", line_t); exit(0);} /* if ((line ==NULL) || (strcmp(line, "@") == 0)) break; */ temp_t = (char *) space(strlen(line_t)+1); (void) sscanf(line_t,"%s",temp_t); int length_t; length_t = (int) strlen(temp_t); for (l = 0; l < length_t; l++) { temp_t[l] = toupper(temp_t[l]); if (!noconv && temp_t[l] == 'T') temp_t[l] = 'U'; } string_t =(char *) space(length_t +11); strcpy(string_t, "NNNNN"); strcat(string_t+5, temp_t); strcat(string_t+5+length_t, "NNNNN"); free(temp_t); char *name_output; name_output=NULL; if(nice){ name_output = (char*) space(sizeof(char)*(length_t+length_s+2)); strcpy(name_output, name_t+1); strcat(name_output, "_"); strcat(name_output, name_s+1); name_output[length_t + length_s +1]='\0'; } if(delta>=0){ snoopT* subopt; snoopT* sub; if(!fast && !plfold_up_flag){ subopt = snoop_subopt(string_t, string_s, delta, 5, penalty, threshloop, threshLE,threshRE, threshDE,threshTE,threshSE,threshD,distance, half_stem, max_half_stem, min_s2, max_s2, min_s1, max_s1, min_d1, min_d2,fullStemEnergy); if(subopt==NULL){printf("no target found under the given constraints\n");free(subopt); free(string_t);free(line_t);continue;} int count=0; for (sub=subopt; sub->structure !=NULL; sub++) { print_struc(sub, string_t, string_s, name_s, name_t, count++,nice); free(sub->structure); } free(subopt); free(string_t); free(line_t); } else if(!plfold_up_flag){ Lsnoop_subopt_list(string_t, string_s, delta, 5, penalty, threshloop, threshLE,threshRE,threshDE,threshTE,threshSE,threshD,distance, half_stem, max_half_stem, min_s2, max_s2, min_s1, max_s1, min_d1, min_d2, alignment_length,name_output,fullStemEnergy); free(string_t); free(line_t); } else { if(plfold_up_flag==1){ int **access_s1; char *file_s1; int s1_len;/* k;*/ /*,j; */ s1_len=strlen(string_t); file_s1 = (char *) space(sizeof(char) * (strlen(name_t+1)+strlen(access)+9)); strcpy(file_s1, access); strcat(file_s1, "/"); strcat(file_s1, name_t+1); strcat(file_s1, "_openen"); access_s1 = read_plfold_i(file_s1, 1, s1_len); if(fast) { Lsnoop_subopt_list_XS(string_t, string_s, (const int**) access_s1, delta, 5, penalty, threshloop, threshLE,threshRE,threshDE,threshTE,threshSE,threshD,distance, half_stem, max_half_stem, min_s2, max_s2, min_s1, max_s1, min_d1, min_d2, alignment_length,name_output, fullStemEnergy); } else{ snoop_subopt_XS(string_t, string_s, (const int **) access_s1, delta, 5, penalty, threshloop, threshLE,threshRE,threshDE,threshTE,threshSE,threshD,distance, half_stem, max_half_stem, min_s2, max_s2, min_s1, max_s1, min_d1, min_d2, alignment_length,name_output,fullStemEnergy); } /* for(j=0; j-1){ free(access_s1[i]); } free(access_s1); } else if(plfold_up_flag==2){ int **access_s1; char *file_s1; int s1_len,k;/* ,j; */ s1_len=strlen(string_t); file_s1 = (char *) space(sizeof(char) * (strlen(name_t+1)+strlen(suffix) + strlen(access)+3)); strcpy(file_s1, access); strcat(file_s1, "/"); strcat(file_s1, name_t+1); strcat(file_s1, "_"); strcat(file_s1, suffix); access_s1 = read_rnaup(file_s1, 1, s1_len); if(fast){ Lsnoop_subopt_list_XS(string_t, string_s, (const int**) access_s1, delta, 5, penalty, threshloop, threshLE,threshRE,threshDE,threshTE,threshSE,threshD,distance, half_stem, max_half_stem, min_s2, max_s2, min_s1, max_s1, min_d1, min_d2, alignment_length, name_output, fullStemEnergy); } else{ snoop_subopt_XS(string_t, string_s, (const int**) access_s1, delta, 5, penalty, threshloop, threshLE,threshRE,threshDE,threshTE,threshSE,threshD,distance, half_stem, max_half_stem, min_s2, max_s2, min_s1, max_s1, min_d1, min_d2, alignment_length,name_output, fullStemEnergy); } /* for(j=0; j-1){ free(access_s1[k]); } free(access_s1); } } } else{ mfe=snoopfold(string_t, string_s, penalty, threshloop, threshLE, threshRE,threshDE, threshD, half_stem, max_half_stem, min_s2, max_s2, min_s1, max_s1, min_d1, min_d2, fullStemEnergy); if(mfe.energy < INF){ print_struc(&mfe, string_t, string_s, name_s, name_t,0,1); free(mfe.structure); } free(line_t); length_t= (int) strlen(string_t); free(string_t); } free(name_t); if(nice){ free(name_output); } } while(1); rewind(mrna); snofree_arrays(strlen(string_s)); /* free's base_pair */ free(string_s);string_s=NULL; free(name_s);name_s=NULL; } while (1); } else{ if(tname==NULL || sname==NULL) RNAsnoop_cmdline_parser_print_help(); sno=fopen(sname, "r"); if(sno==NULL){printf("%s: Wrong snoRNA file name\n", sname);return 0;} mrna=fopen(tname, "r"); if(mrna==NULL){printf("%s: Wrong target file name\n", tname);return 0;} char *temp1[MAX_NUM_NAMES], *temp2[MAX_NUM_NAMES],*AS1[MAX_NUM_NAMES], *AS2[MAX_NUM_NAMES], *names1[MAX_NUM_NAMES], *names2[MAX_NUM_NAMES]; int n_seq, n_seq2; n_seq = read_clustal(mrna, temp1, names1); n_seq2 = read_clustal(sno, temp2, names2); if (n_seq != n_seq2){ for (i=0; temp1[i]; i++) { free(temp1[i]); free(temp2[i]); } nrerror("unequal number of seqs in alignments"); } for(i=0;temp1[i];i++){ AS1[i] = (char*) space((strlen(temp1[i])+11)*sizeof(char)); AS2[i] = (char*) space((strlen(temp2[i])+11)*sizeof(char)); strcpy(AS1[i],"NNNNN"); strcpy(AS2[i],"NNNNN"); strcat(AS1[i],temp1[i]); strcat(AS2[i],temp2[i]); strcat(AS1[i],"NNNNN\0"); strcat(AS2[i],"NNNNN\0"); } for (i=0; temp1[i]; i++) { free(temp1[i]); free(temp2[i]); } AS1[n_seq]=NULL; AS2[n_seq]=NULL; update_fold_params(); alisnofold((const char **)AS2, max_asymm, threshloop, min_s2, max_s2, half_stem, max_half_stem); snoopT struc; struc = alisnoopfold((const char**) AS1, (const char**) AS2, penalty, threshloop, threshLE,threshRE,threshDE,threshD, half_stem, max_half_stem, min_s2, max_s2, min_s1, max_s1, min_d1, min_d2); if(!(struc.structure==NULL)){ aliprint_struc(&struc, (const char**) AS1,(const char**) AS2, names1, names2,0,nice); free(struc.structure); } snoopT* subopt = alisnoop_subopt((const char**) AS1, (const char**) AS2, delta, 5, penalty, threshloop, threshLE,threshRE, threshDE, threshD, threshTE,threshSE,distance, half_stem, max_half_stem, min_s2, max_s2, min_s1, max_s1, min_d1, min_d2); snoopT *sub; if(subopt==NULL){ printf("no target found under the given constraints\n"); }else{ int count=1; for (sub=subopt; !(sub->structure ==NULL); sub++) { aliprint_struc(sub, (const char**) AS1,(const char**) AS2,names1,names2,count,nice); free(sub->structure); count++; } } alisnofree_arrays(strlen(AS2[0])); free(subopt); for (i=0; temp1[i]; i++) { free(AS1[i]); free(AS2[i]); free(names1[i]); free(names2[i]); } } fclose(sno); fclose(mrna); return 0; } static void print_struc(snoopT *dup, const char *s1, const char *s2, char *name_s, char *name_t, int count, int nice) { int l1; l1 = strchr(dup->structure, '&')-dup->structure; char *target_struct; int shift=0, n2; char* psoutput; psoutput = (char*) space(100*sizeof(char)); /* if(dup->i > strlen(s1)-10){ */ /* dup->i--; */ /* l1--; */ /* } */ /* if(dup->i-l1< 0){ */ /* l1--; */ /* shift++; */ /* } */ target_struct = (char*) space(sizeof(char) * (strlen(dup->structure)+1)); strncpy(target_struct, dup->structure+shift, l1); strncat(target_struct, dup->structure + (strchr(dup->structure, '&')-dup->structure), strlen(dup->structure) - (strchr(dup->structure, '&')-dup->structure)); strcat(target_struct,"\0"); char *target; target = (char *) space(l1+1); strncpy(target, (s1+dup->i-l1+5), l1); target[l1]='\0'; char *s4; n2 = strlen(s2); s4 = (char*) space(sizeof(char) *(n2-9)); strncpy(s4, s2+5, n2-10); s4[n2-10]='\0'; printf("%s %3d,%-3d;%3d : %3d,%-3d (%5.2f = %5.2f + %5.2f + %5.2f + %5.2f + 4.1 ) (%5.2f) \n%s&%s\n", target_struct, dup->i+1-l1, dup->i, dup->u, dup->j+1, dup->j + (strrchr(dup->structure,'>') - strchr(dup->structure, '>'))+1, (dup->Loop_D + dup->Duplex_El + dup->Duplex_Er + dup->Loop_E) + 4.10, dup->Duplex_El, dup->Duplex_Er, dup->Loop_E, dup->Loop_D,dup->fullStemEnergy,target,s4); if(nice){ char *temp_seq; char *temp_struc; temp_seq = (char*) space(sizeof(char)*(l1+n2-9)); temp_struc = (char*) space(sizeof(char)*(l1+n2-9)); strcpy(temp_seq, target); strcat(temp_seq, s4); strncpy(temp_struc, target_struct, l1); strcat(temp_struc, target_struct+l1+1); temp_seq[n2+l1-10]='\0'; temp_struc[n2+l1-10]='\0'; cut_point = l1+1; char str[16];char upos[16]; char *temp; int length_name = strlen(name_t)+strlen(name_s); temp = (char *) space(sizeof(char)*(length_name+1)); strcpy(temp,name_t+1); strcat(temp,"_"); strcat(temp,name_s + 1); temp[length_name] = '\0'; strcpy(psoutput,"sno_"); sprintf(str,"%d",count); strcat(psoutput,str); sprintf(upos,"%d",dup->u); strcat(psoutput,"_u_"); strcat(psoutput,upos); strcat(psoutput,"_"); strcat(psoutput,temp); strcat(psoutput,".ps\0"); PS_rna_plot_snoop_a(temp_seq, temp_struc, psoutput, NULL, NULL); cut_point = -1; free(temp_seq);free(temp_struc);free(psoutput);free(temp); } free(s4); free(target_struct); free(target); } static void aliprint_struc(snoopT *dup, const char **s1, const char **s2, char **name_t, char **name_s,int count,int nice) { int l1; l1 = strchr(dup->structure, '&')-dup->structure; char *target_struct; int shift=0; /* if(dup->i > strlen(s1[0])-10){ */ /* dup->i--; */ /* l1--; */ /* } */ /* if(dup->i-l1< 0){ */ /* l1--; */ /* shift++; */ /* } */ int length_struct = strlen(dup->structure); target_struct = (char*) space(sizeof(char) * (length_struct+1)); strncpy(target_struct, dup->structure+shift, l1); strncat(target_struct, dup->structure + (strchr(dup->structure, '&')-dup->structure), length_struct - (strchr(dup->structure, '&')-dup->structure)); strcat(target_struct,"\0"); /* get the corresponding alignment slice */ int n_seq,s; for(s=0; s1[s]!=NULL; s++); n_seq=s; int n1,n2; n1=strlen(s1[0]);n2=strlen(s2[0]); char **target; target = (char**) space((n_seq+1)*sizeof(char*)); for(s=0;si-l1+5), l1); strcat(target[s],"&"); strncat(target[s], (s2[s]+5), n2-10); target[s][l1+n2-9]='\0'; } char * consens; consens = consens_mis((const char**)target); consens[l1]='&'; printf("%s %3d,%-3d;%3d : %3d,%-3d (%5.2f = %5.2f + %5.2f + %5.2f + %5.2f + 4.1; duplex cov = %5.2f; stem cov = %5.2f )\n%s\n", dup->structure, dup->i+1-l1, dup->i, dup->u, dup->j+1, dup->j + (strrchr(dup->structure,'>') - strchr(dup->structure, '>'))+1, (dup->Loop_D + dup->Duplex_El + dup->Duplex_Er + dup->Loop_E)/n_seq + 4.10, dup->Duplex_El/n_seq, dup->Duplex_Er/n_seq, dup->Loop_E/n_seq, dup->Loop_D/n_seq,dup->pscd/n_seq, dup->psct/n_seq,consens ); if(nice){ char* psoutput; psoutput = (char*) space(100*sizeof(char)); char *temp_seq, *temp_struct, **temp_target; temp_seq = (char*) space((l1 + n2 -9)*sizeof(char)); temp_struct = (char*) space((l1 + n2 -9)*sizeof(char)); temp_target = (char**) space((n_seq+1)*sizeof(char*)); for(s=0;si-l1+5), l1); strncat(temp_target[s], (s2[s]+5), n2-10); temp_target[s][l1+n2-10]='\0'; } strncpy(temp_seq, consens, l1); strncpy(temp_struct, target_struct, l1); strcat(temp_seq, consens+l1+1); strcat(temp_struct, target_struct+l1+1); temp_seq[n2-10+l1]='\0'; temp_struct[n2-10+l1]='\0'; char str[16]; char upos[16]; char *temp; int length_name = strlen(name_t[0]) + strlen(name_s[0])+1; temp = (char *) space(sizeof(char)*(length_name+1)); strcpy(temp,name_t[0]); strcat(temp,"_"); strcat(temp, name_s[0]); temp[length_name] = '\0'; strcpy(psoutput,"snoaln_"); sprintf(str,"%d",count); strcat(psoutput,str); sprintf(upos,"%d",dup->u); strcat(psoutput,"_u_"); strcat(psoutput,upos); strcat(psoutput,"_"); for(s=0; sstructure, psoutput, (const char**) target, (const char**) name_t); psoutput[1]='t'; psoutput[2]='r'; PS_rna_plot_snoop_a(temp_seq, temp_struct, psoutput, NULL, (const char**) temp_target); cut_point = -1; free(psoutput); for(s=0; s')){ fprintf(stderr,"file %s is not in RNAup format\n",fname); exit(EXIT_FAILURE); } while(!strstr(fgets(tmp,sizeof(tmp),in), "pos")){}; /* (void) fgets(tmp,sizeof(tmp),in); //Datum */ /* (void) fgets(tmp,sizeof(tmp),in); //white line */ /* (void) fgets(tmp,sizeof(tmp),in); //RNAup */ /* (void) fgets(tmp,sizeof(tmp),in); //sequence length */ /* (void) fgets(tmp,sizeof(tmp),in); //sequence */ int dim_x; dim_x = get_max_u(tmp,'S'); /* get unpaired regions by conting tabs in second line */ access = (int**) space(sizeof(int *) * (dim_x+2)); for(i=0; i< dim_x+2; i++){ access[i] =(int *) space(sizeof(int) * (end_r-beg_r+7)); } for(i=0;i -1) /* read a record */ { float n; /* int i; */ /* int u; */ beg_r--; if(beg_r<1){ if(sscanf(tmp,"%d%n",&seq_pos,&temp)==1){ /* read sequenz pos = 1. spalte */ offset+=temp; int count; count=1; while(sscanf(tmp+offset,"%f%n",&n,&temp)==1){ /* read columns */ offset+=temp; access[count][-beg_r+5+1]=(int) rint(100*n); /* seq_pos+5 */ /* printf("%d %d %f\n", count, -beg_r, access[count][-beg_r]); */ count++; } offset=0; } } } fclose(in); return access; } static int ** read_plfold_i(char *fname, const int beg, const int end) { FILE *in=fopen(fname,"r"); int i,j; int **access; int offset,temp; temp=0;offset=0; int seq_pos; int beg_r, end_r; beg_r=beg; end_r=end; if(in==NULL){ perror(" open error"); exit(EXIT_FAILURE); } char tmp[2048]={0x0}; /* char *ptr; */ if(fgets(tmp,sizeof(tmp),in)==0){ perror("Empty File"); } if(strchr(tmp,'>')){ fprintf(stderr,"file %s is not in RNAplfold format",fname); exit(EXIT_FAILURE); } if(fgets(tmp,sizeof(tmp),in)==0){ perror("No accessibility data"); } int dim_x; dim_x=get_max_u(tmp,'\t'); access = (int**) space(sizeof(int *) * (dim_x+2)); for(i=0; i< dim_x+2; i++){ access[i] =(int *) space(sizeof(int) * (end_r-beg_r+1)); } for(i=0;i -1) /* read a record */ { float n; /* int i; */ /* int u; */ beg_r--; if(beg_r<1){ if(sscanf(tmp,"%d\t%n",&seq_pos,&temp)==1){ offset+=temp; int count; count=1; while(sscanf(tmp+offset,"%f%n",&n,&temp)==1 ){ offset+=temp; access[count][seq_pos+5]=(int) rint( 100 * n); /* round the number */ count++; } offset=0; } } } fclose(in); return access; } static void redraw_output(char *fname, char *output, int plfold_up_flag, char *suffix, char *access, char *sname, char *tname) { char *line; line=NULL; int count; short two_seq = 0; char *results; char *sequence; char *query; char *target; char *structure; char *pos; int posi; int begin=0, end=0, u=0; if(output==NULL){ output = (char*) space(sizeof(char)*2); strcpy(output,".\0"); } count=0; if(sname==NULL && tname==NULL){ while((line=get_line(stdin))!=NULL) { count++; if(two_seq==0 && *line =='>'){ query=(char*) space(strlen(line)+1); (void) sscanf(line,"%s",query); free(line); line=NULL; memmove(query, query+1, strlen(query)); two_seq++; } else if(two_seq==1 && *line =='>'){ target=(char*) space(strlen(line)+1); (void) sscanf(line,"%s",target); free(line); line=NULL; memmove(target, target+1, strlen(target)); two_seq++; } else if(two_seq == 2) { int *relative_access; relative_access=NULL; printf("%s %s\n", target,query); if(strchr(line, '(') && strchr(line,'&') && strchr(line, '(')&& strchr(line, ',') && strchr(line,':') && strchr(line, '-')){ results=line; int length; pos = strchr(results, ' '); posi = (int) (pos - results); length = posi; structure = (char *) space((length+1) * sizeof(char)); sscanf(results,"%s",structure); /* parse structure */ char *line2; if((line2=get_line(stdin))!=NULL){ sequence=(char *) space((length+1)* sizeof(char)); sscanf(line2,"%s",sequence); if(line2!=NULL){ free(line2); } } else{ printf("your file is fucked up");exit(0); } /* parse sequence */ sscanf(pos, "%10d,%10d;%10d", &begin,&end,&u); /* parse coordinates */ if(plfold_up_flag){ int **access_s1; char *file_s1; /* read_rnaup_file */ file_s1 = (char*) space(sizeof(char) * (strlen(target)+strlen(suffix)+strlen(access)+3)); strcpy(file_s1, access); strcat(file_s1, "/"); strcat(file_s1, target); strcat(file_s1, "_"); strcat(file_s1, suffix); access_s1 = read_rnaup(file_s1, begin, end); free(file_s1); relative_access = (int*) space(sizeof(int)*(end-begin+2)); relative_access[0]=access_s1[1][1+5]; int i; for(i=2;i<(end-begin+2);i++){ relative_access[i-1]=access_s1[i+1][i+5] - access_s1[i][i+4]; } int l = access_s1[0][0]; while(--l>-1){ free(access_s1[l]); } free(access_s1); } char *catseq, *catstruct,*output_file; catseq = (char*) space(strlen(sequence) *sizeof(char)); catstruct = (char*) space(strlen(structure) *sizeof(char)); int l1 = strchr(structure, '&')-structure; strncpy(catseq, sequence, l1); strcat(catseq, sequence+l1+1); strncpy(catstruct, structure, l1); strcat(catstruct, structure+l1+1); strcat(catseq,"\0"); strcat(catstruct,"\0"); /* printf("%s\n%s\n%s\n%s", catseq,sequence,catstruct,structure); */ cut_point=l1+1; output_file = (char*) space((strlen(output) + strlen(query)+strlen(target)+50)*sizeof(char)); strcpy(output_file,output); strcat(output_file, "/"); strcat(output_file,"sno_"); strcat(output_file,query); strcat(output_file,"_"); strcat(output_file, target); strcat(output_file,"_u_"); char str[9]; sprintf(str,"%d",u); strcat(output_file,str); strcat(output_file,"_"); sprintf(str,"%d",count); strcat(output_file,str); strcat(output_file,".ps\0"); PS_rna_plot_snoop_a(catseq, catstruct, output_file,relative_access,NULL); if(relative_access){ free(relative_access); } free(output_file);output_file=NULL; free(catseq);catseq=NULL; free(catstruct);catstruct=NULL; free(structure); free(sequence); free(line); line=NULL; } else if(*line=='>'){ free(query);free(target); two_seq=1; query = (char*) space(sizeof(char)*(strlen(line)+1)); (void) sscanf(line,"%s",query); free(line); memmove(query, query+1, strlen(query)); } } } free(target); free(query); } else if(!(tname==NULL && sname==NULL)) { FILE* sno, *mrna; int i; sno=fopen(sname, "r"); if(sno==NULL){printf("%s: Wrong snoRNA file name\n", sname);} mrna=fopen(tname, "r"); if(mrna==NULL){printf("%s: Wrong target file name\n", tname);} char *temp1[MAX_NUM_NAMES], *temp2[MAX_NUM_NAMES], *AS1[MAX_NUM_NAMES], *AS2[MAX_NUM_NAMES], *names1[MAX_NUM_NAMES], *names2[MAX_NUM_NAMES]; int n_seq, n_seq2; n_seq = read_clustal(mrna, temp1, names1); n_seq2 = read_clustal(sno, temp2, names2); if (n_seq != n_seq2){ for (i=0; temp1[i]; i++) { free(temp1[i]); free(temp2[i]); } nrerror("unequal number of seqs in alignments"); } for(i=0;temp1[i];i++){ AS1[i] = (char*) space((strlen(temp1[i])+11)*sizeof(char)); AS2[i] = (char*) space((strlen(temp2[i])+11)*sizeof(char)); strcpy(AS1[i],"NNNNN"); strcpy(AS2[i],"NNNNN"); strcat(AS1[i],temp1[i]); strcat(AS2[i],temp2[i]); strcat(AS1[i],"NNNNN\0"); strcat(AS2[i],"NNNNN\0"); } for (i=0; temp1[i]; i++) { free(temp1[i]); free(temp2[i]); } AS1[n_seq]=NULL; AS2[n_seq]=NULL; int count=0; while((line=get_line(stdin))!=NULL) { results=line; if(strchr(line, '(') && strchr(line,'&') && strchr(line, '(')&& strchr(line, ',') && strchr(line,':') && strchr(line, '-')) { count++; int length; /* int sbegin, send; */ /* int energy; */ snoopT dup; pos = strchr(results, ' '); posi = (int) (pos - results); length = posi; dup.structure = (char *) space((length+1) * sizeof(char)); sscanf(results,"%s %10d,%10d;%10d : %10d,%10d (%10f = %10f + %10f + %10f + %10f + 4.1; duplex cov = %10f; stem cov = %10f", dup.structure, &begin, &(dup.i), &(dup.u), &(dup.j), &end, &(dup.energy), &(dup.Duplex_El), &(dup.Duplex_Er), &(dup.Loop_E), &(dup.Loop_D), &(dup.pscd), &(dup.psct)); /* parse duplex stuff; */ dup.energy*=n_seq; dup.Duplex_El*=n_seq; dup.Duplex_Er*=n_seq; dup.Loop_E*=n_seq; dup.Loop_D*=n_seq; aliprint_struc(&dup, (const char**) AS1, (const char**)AS2, names1, names2,count,1); free(dup.structure); free(line); } else{ free(line); } } for (i=0; AS1[i]; i++) { free(AS1[i]); free(AS2[i]); free(names1[i]); free(names2[i]); } fclose(mrna); fclose(sno); } /* free(output); */ }