/* Code by Ethan Wolf 1995. Graphics code by David Wilson. */ /*** See dimer_trimer_like.c, function calculate_likelihood to see how */ /*** to compute likelihood line number before call DimerTrimerScore() */ #define XCOORD_SCORED_BY 430 #define XCOORD_SEQSCORES 0 #define YCOORD_SEQSCORES 55 /* Each char takes "15" so put us below printing */ #define WIDTH_ADD_ON 30 /** Empty space after seq in main window. **/ #define xmin 120 /* Want this far enough over so have room for seqscoreb */ /** seqscoreb is 108 wide. xmin used to be 60 **/ #include "sc.h" #include #include "sc2seq.h" /* This also includes scconst.h, sc.h, interface.h */ #include #include #include #include "scio.h" #include "options.h" Display *display; Window window, zoom, printb, quit, next, prev, autob, zoomb, methodb, paramb, helpb, help, seqscoreb; GC gc, gcgray, gcdash; #include #include #include #include #include "config.h" #include "switches.h" #include "hmm.h" /* FILE *debug_out; */ /*************************** MULTICOIL LIKELIHOOD STUFF. *********************/ double both_determinants[2][MAX_CLASS_NUMBER]; double *determinants= &both_determinants[0][0]; double both_inverse_covars[2][MAX_CLASS_NUMBER][MAX_NUM_SCORE_DIM][MAX_NUM_SCORE_DIM]; double *inverse_covars= &both_inverse_covars[0][0][0][0]; double both_covars_submatrix[2][MAX_CLASS_NUMBER][MAX_NUM_SCORE_DIM][MAX_NUM_SCORE_DIM]; double both_means_submatrix[2][MAX_CLASS_NUMBER][MAX_NUM_SCORE_DIM]; double *means_submatrix= &both_means_submatrix[0][0][0]; double init_class_prob[MAX_CLASS_NUMBER]; int multi_functnum[MAX_TABLE_NUMBER]; int number_score_dim=0; int number_classes=0; int number_multi_lib[MAX_TABLE_NUMBER]; /* Used to detrmine the number of different library */ /* distance values to do gaussian param over. */ #ifdef TEST_VERSION FILE *ftotal_like=NULL; /* Extern for now. ***/ double total_gauss_like[GRID_SIZE_DIMER][GRID_SIZE_TRIMER]; /* When ftotal_like is read from get_defaults(), the */ /* i,j entry will hold the total gaussian like when */ /* assuming dimers have init_prob = */ /* MIN_DIMER_PROB + i*DIMER_GRID_STEP */ /* and trimers have init_prob = */ /* MIN_TRIMER_PROB + j*TRIMER_GRID_STEP */ /* These constants are defined in scconst.h. */ /* The value is updated during "output_seq()". */ /*************************** PAIRDIFF LIKELIHOOD STUFF. *********************/ double differ_determinants[MAX_CLASS_NUMBER]; double differ_inverse_covars[MAX_CLASS_NUMBER][MAX_NUM_SCORE_DIM][MAX_NUM_SCORE_DIM]; double differ_covars_submatrix[MAX_CLASS_NUMBER][MAX_NUM_SCORE_DIM][MAX_NUM_SCORE_DIM]; double differ_means_submatrix[MAX_CLASS_NUMBER][MAX_NUM_SCORE_DIM]; int number_differ_score_dim =0; int number_differ_lib=0; int differ_window_length = 28; int which_differentiator=1; /** Default is just the score for the window **/ /** 2 means to give each residue a pair of scores**/ /** corresponding to the max and min of the **/ /** windows containing it. **/ int differ_functnum=0; #endif /************************** General Scoring parameters. ********************/ int functnum=0, pair_functnum=0; int preproc_table=0; double log_bound = 0; /*In PRN_MODE only scores > bound output to logfile.*/ int WINDOW = PAIRWINDOW; int window_length[MAX_TABLE_NUMBER]; double SCALE0 = DEFAULT_SCALE0; double scale0s[MAX_TABLE_NUMBER]; /** Both used as bonnie's scale0 and **/ double scale0p[MAX_TABLE_NUMBER]; /** in computing posterior prob. from **/ /** prior as sigma^2(prior)=scale0*sigma^2 */ /** where these are the standard deviation * /** of the n*(prior_prob) and the distib */ /** of number of hits when sampling */ int log_offset_to_use; int offset_to_use=-1; /* Default offset in DimerTrimerScore is all. */ int avg_max=1; /* Default means use max coil scores in log file. */ /** Value of 0 means do average, 1 means do max. **/ /** 2 means do both. */ int number_tables =0; int method = MultiCoil; int main_method = MultiCoil; int main_preprocessor_method = MultiCoil; int preprocessor_method = MultiCoil; int main_table = 0; int table = 0; /*********************** Parameters involved in Coil Scores ***********/ int by_coil_or_seq =0; /*** To determe what type of coil_out file to do. **/ int weighted_avg = 0; /* Determines which method for determining coil */ /* score will be used in DimerTrimerScore. */ int start_coil_at_reg_shift =0; /* For averaging scores over coil length. */ int Coil_Score =0; /* 0 means for MultiCoil do residues score, 1 means */ /* average scores over coils, 2 means take max coil */ /* score. */ Sequence sequence; char *methodname[] = {"MultiCoil", "Coil","NEWCOILS", "HMMCoil","ActualCoil", "PairCoilDiff", "PairCoilDiffAvg"}; /* "Coil" is PairCoil or */ char table_diff_stat_filename[MAXLINE]="/dev/null"; char weightfile_outname[MAXLINE]="/dev/null"; /***/ void *hmm; /***/ double seq_scores[MAX_CLASS_NUMBER]; double *all_scores; double *scores; /** The space these pointers points to is alocated by **/ /** using static variable arrays in ScoreSeq */ /** all_scores is scores[MAX_NUM_SCORE_DIM][MAXSEQLEN][POSNUM +1]*/ /** scores is scores[current_score_algor][MAXSEQLEN][POSNUM+1] **/ double all_preprocess_like[MAX_NUM_SCORE_DIM][MAXSEQLEN][POSNUM+1]; /** Since preprocess_likelihood is used everywhere, just */ /** make it an external variable. */ /** This is given a value in ScoreSeq and sometime in */ /** ActualCoil(). */ double *preprocess_like; /* This will point to the current table version */ /* of all_preprocess_like */ char all_offsets[MAX_NUM_SCORE_DIM][MAXSEQLEN]; /************************ XWindow stuff. ***********************************/ int zoomptr=0; int ps_res_per_line=0; /* 0 is a flag that it wasn't set in get_defaults */ /* so should just print it using nominal_length */ /* as done on the screen. If it is set to -1 */ /* then should print it all on one line. */ char *motd[] = { " *********************************", " ** **", " ** MultiCoil (version 1.0) **", " ** **", " ** **", " *********************************", "", "MultiCoil is described in:", " Ethan Wolf, Peter Kim, Bonnie Berger,", " `MultiCoil: A Program for Predicting Two- and " Three-stranded Coiled Coils'.", " ", "Email inquiries to: multicoil-help@theory.lcs.mit.edu.", "" }; char *hmotd[] = { "help: Clicking on help displays this message.", " ", "next: Displays the the next sequence.", " Same as pressing 'n' or down-arrow.", " ", "prev: Displays the previous sequence.", " Same as pressing 'p' or up-arrow.", " ", "auto: Automatically steps through the sequences", " of a file.", " ", "zoom: Gives more information on the highest", " scoring residue. Same as pressing 'z'.", " Clicking on a specific region will give", " more information on that region.", " ", "print: Prints the sequence.", " ", /********** "table: Switches the table used by PairCoil and", " by TableDiff methods, same as 't'. Note that", " key 'i' (input_table) changes the preprocessor's table.", " ", **********/ "method: Switches the scoring method, same as 'm'.", " 'w' swithces Which preprocessor is used. ", " ", "key 'Up' or 'Down' moves the bound up or down by .01", "key 'Return' rescores the sequence after changing the bound", " ", "key 's' or 'l': switch between singles and pair", " probabilities at lib distances.", " ", "key 'o': Calc paircoil class prob. with respect only to", " positive oligimerization classes (ignore pdb- class)", "key '7': MultiCoil max score is obtained by taking tuple", " of each dimension's max score (combo register).", "key 'x': MultiCoil computes scores for each offset, taking", " max of those scores to be the max score (max reg).", "key 'r': Cycle through offsets 0-7 for register (7=combo reg)", "", "key 'c': Switch between MultiCoil register score, average coil", " score, and max coil score", "key 'b': For Coil score, toggle whether start new coils at", " a register shift (break in structure).", "key 'a': For Coil score, toggle between unweighted and weighted avg", "", "quit: Quits the program and completes log file.", " 'q' quits without completing the log.", " Pressing 'q' in any window also closes ", " that window.", "" }; ShowTitle (Window win) { int i; for (i = 0; motd[i][0]; ++i) XDrawString(display, win, gc, 5, 17*i+20, motd[i], strlen(motd[i])); XDrawString(display,win,gcgray,100,17*i+20, "Press Any Key To Continue",25); } ShowHelp (Window win) { int i; /* Places hmotd text in win */ for (i=0; hmotd[i][0]; ++i) XDrawString(display,win,gc,5,17*i+20, hmotd[i],strlen(hmotd[i])); } static int seqnum=0, sequp=0, seqlow=1; char autoadvance=0; main (int argc, char *argv[]) { struct itimerval interval; fd_set fds, fdnone; FILE *fin, *ps; XEvent report; char buff[MAXLINE+MAXSEQLEN]; /* Note that fout_raw[MAX_TABLE_NUMBER+1] and fout_coils[MAX_TABLE_NUMBER+1] */ /* will be files of tuples of each residue's scores according to all tables. */ /* Never mind... **/ char *command_line = NULL; FILE *fgin=NULL, *fout=NULL, *flog=NULL, *flog_coil_conflicts=NULL, *fout_coils=NULL; FILE *fpin[MAX_TABLE_NUMBER]; FILE *pir=NULL; char pir_name[MAXLINE]; char likelihoods[MAX_TABLE_NUMBER][MAXLINE]; /* Files of likelihood lines. */ char printfile[500]; /* The filename to output postscript print to. */ char *print= printfile; extern int number_classes; char class_sc_filenames[2][MAX_CLASS_NUMBER][MAXLINE]; double both_class_covars[2][MAX_CLASS_NUMBER][MAX_NUM_SCORE_DIM] double both_class_means[2][MAX_CLASS_NUMBER][MAX_NUM_SCORE_DIM]; [MAX_NUM_SCORE_DIM]; char gauss_param[2][MAXLINE]; /** 0 entry is for combo offset, 1 entry is for **/ /** max of other registers. ***/ extern int functnum; /* What is used in scscore.c. Takes value from */ /* pair_functnum and multi_functnum. */ char lib[MAXFUNCTNUM]; /* Use MAXFUNCTUM=7 here since are 7 distances */ /* in coil, so we know AT MOST 7 lib functions. */ extern int pair_functnum; /* The number of distances in lib[].*/ char multi_lib[MAX_TABLE_NUMBER][MAXFUNCTNUM]; extern int multi_functnum[MAX_TABLE_NUMBER]; #ifdef TEST_VERSION char differ_class_sc_filenames[MAX_CLASS_NUMBER][MAXLINE]; double differ_means[MAX_CLASS_NUMBER][MAX_NUM_SCORE_DIM]; double differ_covars[MAX_CLASS_NUMBER][MAX_NUM_SCORE_DIM][MAX_NUM_SCORE_DIM]; char differ_gauss_param[MAXLINE]; char differ_lib[MAXFUNCTNUM]; /* Distance values for PairCoilDiff */ extern int differ_functnum; #ifdef TEST_VERSION int acnt = 1; int mode = 0; /* If the i-th option is set, mode gets 1 in i-th bit. */ double bound =0; double upper_bound; /* In DEBUG_MODE (for printing out the predicted */ /* heptad repeat), scores in [bound,upper_bound] */ /* are printed in capital letters. */ double m[MAX_TABLE_NUMBER], b[MAX_TABLE_NUMBER]; /*The likelihood line is mx+b for paircoil. */ double m_single[MAX_TABLE_NUMBER], b_single[MAX_TABLE_NUMBER]; /* Like line for singlecoil. */ int i, class; int table_to_remove_from = -1; /* Used if want to remove ALL input seq */ /* from that table at outset. **/ int old_num_dim_table[2]; int new_num_dim_table[2]; int old_num_dim_differ =0; /* The gauss param file for multiple distances is */ /* contains the parameters for all distances. */ #ifdef TEST_VERSION double prior_freq_single[MAX_TABLE_NUMBER], prior_freq_pair[MAX_TABLE_NUMBER]; int which_priors =0, which_priorp =0; int good_turing_fixed_disc =0; /** If set to 1 do the good_turing_fixed_disc estimate of probs*/ int structural_pos[POSNUM+1]; /** Holds list of which positions are **/ /** structurally important. All others **/ /** get genbnk probabilities.... **/ /** Value -1 indicates no more entries. **/ #endif /*** INITIALIZATION ***/ fprintf(stderr,"\nAbout to initialize"); initialize(&number_classes, &number_multi_lib[0], window_length, scale0s, scale0p,fpin[0], &likelihoods[0][0], &pir_name[0], print, gauss_param, differ_gauss_param, &mode, prior_freq_single, prior_freq_pair, structural_pos); fprintf(stderr,"\nInitialized"); read_command_line(argc, argv, &command_line, multi_functnum, &mode, number_multi_lib, multi_lib); fprintf(stderr,"\n\nRead commandline, command_line=%s",command_line); /* Get the default options in the .paircoil file!!! */ get_defaults(argv[1], &mode, &log_bound, &fgin, fpin, &ftotal_like, &fout, &flog, &flog_coil_conflicts, &fout_coils, &by_coil_or_seq, likelihoods, pir_name, lib, differ_lib, multi_lib, &pair_functnum, &differ_functnum, multi_functnum, print, &main_method, &main_preprocessor_method, &main_table, &number_tables, &offset_to_use, &avg_max, &weighted_avg, &Coil_Score, &start_coil_at_reg_shift, table_diff_stat_filename, weightfile_outname, &ps_res_per_line, &number_classes, class_sc_filenames, differ_class_sc_filenames, gauss_param, differ_gauss_param, number_multi_lib, &number_differ_lib, init_class_prob, &table_to_remove_from, command_line, window_length, scale0s, scale0p, prior_freq_single, prior_freq_pair, &which_priors, &which_priorp, &good_turing_fixed_disc, structural_pos, &which_differentiator, &differ_window_length); fprintf(stderr,"\nOut of get_def"); bound = log_bound; #ifdef TEST_VERSION compute_multivariate_num_differ_dim(&number_differ_score_dim, number_differ_lib, differ_functnum, &old_num_dim_differ, which_differentiator); #endif compute_multivariate_num_dim(mode, &number_score_dim, number_multi_lib,multi_functnum, number_tables, new_num_dim_table, old_num_dim_table); if (number_classes >0) { /*** Get gaussian param for the classes. */ if (ftotal_like) { /** Want to output total gauss like over classes*/ for (i=0; i1?argv[1]:"-","r"); fprintf(stderr, "\nfin is opened to %s",argc>1?argv[1]:"-"); /*** Note: If can't open xdisplay then OpenScreen will write log and */ /*** out files if those options are set. */ fprintf(stderr,"\nInitializing\n"); OpenScreen(mode,log_bound,fin,flog,flog_coil_conflicts, fout_coils, fout, m,b,m_single, b_single, lib,differ_lib,multi_lib, main_method,main_preprocessor_method, &seqnum, which_priors, which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); painter(); /* Let user read message while initializing */ FD_ZERO(&fds); FD_ZERO(&fdnone); interval.it_value.tv_sec = 2; interval.it_value.tv_usec = 0; interval.it_interval.tv_sec = 2; interval.it_interval.tv_usec = 0; nopainter(); XSelectInput(display, window, KeyPressMask | ButtonPressMask | ExposureMask); while (1) { if (autoadvance) do { FD_SET(ConnectionNumber(display),&fds); if (!select(FD_SETSIZE,&fds,&fdnone,&fdnone,&interval)) NextSeq(fin,lib,differ_lib,multi_lib, m,b,m_single,b_single, mode,bound,flog,flog_coil_conflicts,fout_coils, fout, which_priors, which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); } while ( (!XPending(display)) && (autoadvance)); /* Autoadvance until button is hit or end of seq file. */ XNextEvent(display,&report); switch (report.type) { case Expose: if (report.xany.window == zoom) Zoomer(); else if (report.xany.window == help) ShowHelp(help); else if (sequence.seqlen) ShowScores(scores,&sequence,NULL,bound, mode); else ShowTitle(window); break; case ButtonPress: /* if (!sequence.seqlen) {NextSeq(fin,lib,differ_lib,multi_lib, m,b,m_single, b_single, bound,flog,flog_coil_conflicts, fout_coils,fout, which_priors, which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); break;} */ if (report.xany.window == help) {XUnmapWindow(display,help); break;} if (report.xany.window == zoom) {XUnmapWindow (display, zoom); break;} if (report.xany.window == window) { /* Zoom in and show register */ int tmp; tmp = XYtoPos(report.xbutton.x,report.xbutton.y); if (tmp) { ZoomPnt(); zoomptr=tmp; Zoomer(); ZoomPnt(); } } if (report.xany.window == helpb) XMapWindow(display,help); if (report.xany.window == next) NextSeq(fin,lib,differ_lib,multi_lib, m,b,m_single, b_single, mode,bound,flog,flog_coil_conflicts,fout_coils,fout, which_priors,which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); if (report.xany.window == prev) PrevSeq(fin,lib,differ_lib,multi_lib, m,b,m_single,b_single, mode,bound, which_priors,which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); if (report.xany.window == autob) /* toggle the automatic sequencer */ { if (autoadvance ^= 1) NextSeq(fin,lib,differ_lib,multi_lib, m,b,m_single,b_single, mode,bound,flog, flog_coil_conflicts,fout_coils, fout, which_priors,which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); } else autoadvance = 0; if (report.xany.window == zoomb) { if (zoomptr) ZoomPnt(); zoomptr = GoodZoom(); Zoomer(); ZoomPnt(); } if (report.xany.window == printb) { XBell(display,-10); ps = sopen(print?print:defprint,"a"); PSDefs(ps, sequence.seqlen); ShowScores(scores,&sequence,ps,bound,mode); sclose(ps); } /*** if (report.xany.window == tableb) { table = (table +1) % number_tables; ShowSeq(lib,differ_lib,multi_lib, m,b,m_single,b_single,mode,0,0, bound); if (zoomptr) Zoomer(); } *****/ if (report.xany.window == methodb) { method_cycle(&method,NUMBER_METHODS); ShowSeq(lib,differ_lib,multi_lib, m,b,m_single,b_single,mode,1,0, bound); if (zoomptr) Zoomer(); } if (report.xany.window == quit) { if (flog || flog_coil_conflicts || fout || ftotal_like || fout_coils) { seqnum=sequp; finish_log(mode,log_bound,fin,flog,flog_coil_conflicts,fout_coils, fout,m,b,m_single,b_single, lib,differ_lib,multi_lib,main_method,main_table, main_preprocessor_method, &seqnum, avg_max, which_priors,which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); } sclose(fin); if (flog) sclose(flog); if (fout) sclose(fout); if (fout_coils) sclose(fout_coils); if (flog_coil_conflicts) sclose(flog_coil_conflicts); if (ftotal_like) sclose(ftotal_like); exit(0); } break; case KeyPress: autoadvance = 0; /* if (!sequence.seqlen) {NextSeq(fin,lib,differ_lib,multi_lib, m,b,m_single, b_single, mode,bound,flog,flog_coil_conflicts, fout_coils,fout, which_priors, which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); break;} */ if (report.xany.window == help) {XUnmapWindow(display,help); break;} switch (XLookupKeysym(&(report.xkey),0)) { case XK_Up: bound = bound + .01; if (bound > 1) bound =1; draw_bound_box(bound, paramb); break; case XK_Down: bound = bound - .01; if (bound < 0) bound =0; draw_bound_box(bound, paramb); break; case XK_Return: /** Signal that bound has been set, should rescore **/ /** -1 then means just need to recompute MultiCoil coil and seq scores. */ ShowSeq(lib,differ_lib,multi_lib, m,b,m_single,b_single,mode,-1,-1, bound); if (zoomptr) Zoomer(); break; case XK_C: /** Switch between register, avg, and max **/ case XK_c: /** coil score for MultiCoil. **/ Coil_Score = (Coil_Score +1) % 3; if (Coil_Score == 0) ShowSeq(lib,differ_lib,multi_lib, /** Change to res score **/ m,b,m_single,b_single,mode,1,1, -2,-2, bound); else ShowSeq(lib,differ_lib,multi_lib, /** Change to coil score with */ m,b,m_single,b_single,mode,-1,-1, bound); /* out recomputing */ if (zoomptr) Zoomer(); /* res score. */ break; case XK_O: /** Oligomerization ratio. ***/ case XK_o: mode ^= ONLY_COILED_CLASSES; /* To XOR to complement that bit. */ if (number_classes >0) { if (method == MultiCoil) type_class_score_convert(sequence,all_scores,bound, mode & ONLY_COILED_CLASSES,number_tables, number_classes); if (preprocessor_method == MultiCoil) type_class_score_convert(sequence,all_preprocess_like,bound, mode & ONLY_COILED_CLASSES,number_tables, number_classes); } ShowSeq(lib,differ_lib,multi_lib, m,b,m_single,b_single,mode,0,2, bound); if (zoomptr) Zoomer(); break; case XK_Print: XBell(display,-10); ps = sopen(print?print:defprint,"w"); PSDefs(ps, sequence.seqlen); ShowScores(scores,&sequence,ps,bound,mode); sclose(ps); XBell(display,-10); break; case XK_Left: if (zoomptr>1) {ZoomPnt(); --zoomptr; Zoomer(); ZoomPnt();} break; case XK_Right: if (zoomptr && zoomptr painter () { struct itimerval value; signal(SIGALRM,Alarm); getitimer(ITIMER_REAL,&value); value.it_value.tv_sec = 0; value.it_value.tv_usec = 500; value.it_interval.tv_sec = 0; value.it_interval.tv_usec = 500; setitimer(ITIMER_REAL,&value,NULL); Alarm(); } nopainter () { struct itimerval value; getitimer(ITIMER_REAL,&value); value.it_value.tv_sec = 0; value.it_value.tv_usec = 0; setitimer(ITIMER_REAL,&value,NULL); XDrawImageString(display,window,gcgray,500,17," ",14); XDrawImageString(display,window,gcgray,225,729,"press any key to continue",25); } ZoomPnt () { int x,y; XPoint pts[3]; if (!zoomptr) return; PostoXY(zoomptr,&x,&y); pts[0].x = x; pts[0].y = y; pts[1].x = -5; pts[1].y = 5; pts[2].x = 10; pts[2].y = 0; XFillPolygon(display,window,gcdash,pts,3,Convex,CoordModePrevious); XDrawLine(display,window,gcdash,x,y+5,x,y+8); } /* Find most likely residue, break ties arbitrarily */ GoodZoom () { int i, zoom = -1; double like=-1; for (i=0; i like) { like = scores[8*i+7]; zoom = i; } return zoom+1; } Zoomer () { int i, offset; double *sc; char buf[32]; XMapWindow(display, zoom); /* Indicate where in the sequence we are */ sprintf(buf,"Position %d ",zoomptr); XDrawImageString(display,zoom,gc,7,20, buf,strlen(buf)); for (i = -3; i<=3; ++i) { if (zoomptr+i>0 && zoomptr+i<=sequence.seqlen) buf[0] = numaa(sequence.seq[zoomptr-1+i]); else buf[0] = ' '; XDrawImageString(display,zoom,gc,159+9*i,20, buf,1); } XDrawLine(display,zoom,gc,159,22,168,22); /* Indicate the probability of being in any given state */ sc = &(scores[8*(zoomptr-1)]); for (i=0; i<=7; ++i) if (sc[i] == sc[7]) {offset=i; break;} /* Used to print out maxscoring */ /* register as 'a' + offset. */ /* Offset 7 will be a flag that no particular offset scored max. */ /* Give likelihood */ sprintf(buf,"%5.1f%% chance",100*sc[7]); XDrawImageString(display,zoom,gc,220,20, buf,strlen(buf)); /***** fprintf(stderr,"\noffset= %d, reg = %c",all_offsets[table][zoomptr-1], 'a' + (zoomptr-1 + all_offsets[table][zoomptr-1])%POSNUM); *******/ if ( sc[7] != 0) if (offset != 7) sprintf(buf,"register %c",'a'+offset); else sprintf(buf,"combo reg "); else sprintf(buf," "); XDrawImageString(display,zoom,gc,230,40, buf,strlen(buf)); for (i=0; i<7; ++i) { sprintf(buf,"%c:%4.1f%% ",'a'+i,(100*sc[i])); XDrawImageString(display,zoom,gc,15+92*(i&3),60+20*(i/4),buf,strlen(buf)); } XDrawString(display,zoom,gc,10,110,"Click on viewgram or use arrow",30); XDrawString(display,zoom,gc,10,127,"keys to view other positions.",29); XDrawString(display,zoom,gc,10,144,"You may need to move this",25); XDrawString(display,zoom,gc,10,161,"window to see the viewgram.",27); } /*********The following is for storing sequences. ***/ static char seqbuff[40000]; static Sequence seqs[30]; static int seqbufpos[30]; static int seqsizes[30]; ReadSeq() { int size, pos, limsize; /* size = seqsizes[sequp-seqnum]; pos = seqbufpos[sequp-seqnum]; */ sequence = seqs[sequp-seqnum]; /** Retrieve seq from memory. **/ /* limsize = (pos+size<=40000) ? size : 40000-pos; bcopy(&(seqbuff[pos]),sequence.code,limsize); bcopy(seqbuff,&(sequence.code[limsize]),size-limsize); */ } StoreSeq() { int size, pos, limsize,i; char *temp, *current; bcopy(seqs,&(seqs[1]),29*sizeof(Sequence)); bcopy(seqbufpos,&(seqbufpos[1]),29*sizeof(int)); bcopy(seqsizes,&(seqsizes[1]),29*sizeof(int)); if (sequp-seqlow>=30) seqlow = sequp-29; /* The 2 accounts for the end of string characters in .code and .title. */ size = 2 + strlen(sequence.code) + strlen(sequence.title) + sequence.seqlen * (sequence.reg ? 2 : 1); pos = seqbufpos[1] + seqsizes[1]; /* Number of characters in array */ /* if (pos>=40000) pos -= 40000;*/ /* limsize = (pos+size<=40000) ? size : 40000-pos; */ if (pos+size >=40000) pos=0; /* Array will be full so start at beginning. */ limsize = (pos+size<=40000) ? size : 40000-pos; /* Save the "code+title+seq+reg" in something permanant, and wrap around at end of array. */ /* David uses bcopy because wants to copy past end of string characters */ /* bcopy(sequence.code,&(seqbuff[pos]),limsize); bcopy(&(sequence.code[limsize]),seqbuff,size-limsize); */ sequence.code=strcpy(&(seqbuff[pos]),sequence.code); sequence.title= strcpy(&(seqbuff[pos+strlen(sequence.code)+1]),sequence.title); temp=sequence.seq; sequence.seq= &(seqbuff[pos+strlen(sequence.code) + strlen(sequence.title) + 2]); for(i=0; i= pos && seqbufpos[sequp-seqlow] < pos+limsize || seqbufpos[sequp-seqlow] < size-limsize)) ++seqlow; seqs[0] = sequence; seqs[0].seqlen= sequence.seqlen; seqsizes[0] = size; seqbufpos[0] = pos; } method_cycle(int *Method, int Number_of_Methods) { do { *Method = (*Method+1)%Number_of_Methods; } while ((*Method == MultiCoil && MultiCoil_Method == unavailable) || (*Method == PairCoil && PairCoil_Method == unavailable) || (*Method == NEWCOILS && NEWCOILS_Method == unavailable) || (*Method == HMMCoil && HMMCoil_Method == unavailable) || /**** (*Method == TableDiffRes && TableDiffRes_Method == unavailable) || (*Method == TableDiffCoil && TableDiffCoil_Method == unavailable) || ****/ (*Method == ActualCoil && ActualCoil_Method == unavailable) || (*Method == PairCoilDiff && PairCoilDiff_Method == unavailable)|| (*Method == PairCoilDiffAvg && PairCoilDiffAvg_Method == unavailable) || (*Method == ActualCoil) && (!sequence.reg) /* Don't do actual coil */ ); /* if didn't input register */ } PrevSeq (FILE *fin, char lib[MAXFUNCTNUM], char differ_lib[MAXFUNCTNUM], char multi_lib[MAX_TABLE_NUMBER][MAXFUNCTNUM], double *m, double *b, double *m_single, double *b_single, int mode, double bound, int which_priors, int which_priorp, double prior_freq_single[MAX_TABLE_NUMBER], double prior_freq_pair[MAX_TABLE_NUMBER], int good_turing_fixed_disc, int structural_pos[POSNUM+1]) { if (seqnum<=seqlow) {XBell(display,0); return;} --seqnum; ReadSeq(); NewSeq(mode,sequence, which_priors, which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); ShowSeq(lib, differ_lib, multi_lib, m, b, m_single, b_single, mode,1,1, bound); } NextSeq (FILE *fin,char lib[MAXFUNCTNUM], char differ_lib[MAXFUNCTNUM], char multi_lib[MAX_TABLE_NUMBER][MAXFUNCTNUM], double *m, double *b, double *m_single, double *b_single, int mode, double bound, FILE *flog, FILE* flog_coil_conflicts, FILE *fout_coils, FILE *fout, int which_priors, int which_priorp, double prior_freq_single[MAX_TABLE_NUMBER], double prior_freq_pair[MAX_TABLE_NUMBER], int good_turing_fixed_disc, int structural_pos[POSNUM+1]) { char title[MAXLINE],code[MAXLINE]; char seq[MAXSEQLEN], reg[MAXSEQLEN]; Sequence temp_seq; int first_time_seq=0; /******** fprintf(stderr,"\nIn next seq"); *********/ if (seqnum < sequp) { ++seqnum; ReadSeq(); } else { temp_seq= sequence; sequence.seqlen = 0; sequence.seq = seq; sequence.reg = reg; sequence.title = title; sequence.code = code; /* if ( ( (mode & POS_MODE) ? !getpos(fin,&sequence) : !getseq2(fin,&sequence) ) ) */ /*** getseq2 now gets register info too... ***/ if (!getseq2(fin,&sequence)) {XBell(display,0); sequence= temp_seq; autoadvance=0; return;} ++seqnum; ++sequp; StoreSeq(); first_time_seq =1; } NewSeq(mode,sequence, which_priors, which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); if (first_time_seq) { output_seq(lib,differ_lib, multi_lib, m,b,m_single, b_single, mode,log_bound,flog,flog_coil_conflicts,fout_coils,fout, avg_max, main_method, main_preprocessor_method, main_table); } ShowSeq(lib, differ_lib,multi_lib, m, b, m_single, b_single,mode,1,1, bound); } /** This is stuff need to always do before rescore a sequence. **/ NewSeq_nonX (int mode, Sequence sequence, int which_priors, int which_priorp, double prior_freq_single[MAX_TABLE_NUMBER], double prior_freq_pair[MAX_TABLE_NUMBER], int good_turing_fixed_disc, int structural_pos[POSNUM+1]) { if (mode & TST_MODE0) recalc_prob(sequence, 0, mode, which_priors, which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); if (mode & TST_MODE1) recalc_prob(sequence,1, mode, which_priors, which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); } NewSeq (int mode, Sequence sequence, int which_priors, int which_priorp, double prior_freq_single[MAX_TABLE_NUMBER], double prior_freq_pair[MAX_TABLE_NUMBER], int good_turing_fixed_disc, int structural_pos[POSNUM+1]) { NewSeq_nonX (mode, sequence, which_priors, which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); XUnmapWindow(display, zoom); zoomptr=0; } /* The likelihood line is mx +b */ void output_seq(char lib[MAXFUNCTNUM],char differ_lib[MAXFUNCTNUM], char multi_lib[MAX_TABLE_NUMBER][MAXFUNCTNUM], double *m,double *b, double *m_single, double *b_single, int mode, double log_bound,FILE *flog, FILE *flog_coil_conflicts, FILE *fout_coils, FILE *fout, int avg_max, int main_method, int main_preprocessor_method, int main_table) { double maxother=0; double maxscore; extern double all_preprocess_like[MAX_NUM_SCORE_DIM][MAXSEQLEN][POSNUM+1]; int i,j; int do_preproc; double junk_scores[MAXSEQLEN][POSNUM+1]; /*** fprintf(stderr,"\nIn outputseq"); ****/ if (is_not_differentiator(main_method)) do_preproc = 0; else do_preproc =1; all_scores = ScoreSeq (lib, differ_lib, multi_lib, m, b, m_single, b_single, mode, main_table, main_method, main_preprocessor_method, log_offset_to_use, &maxscore, 1, do_preproc, log_bound); if (ftotal_like) { add_to_total_likes(all_scores, total_gauss_like, sequence.seqlen, ftotal_like); } /****** fprintf(stderr,"\nGot all_scores"); *******/ scores = &all_scores[main_table* MAXSEQLEN*(POSNUM+1) ]; /* This is all_scores[main_table + 0*number_tables][0][0] */ /* so is score on main_table using library # 0 **/ /* Similarly for preproces_like. **/ preprocess_like = &all_preprocess_like[main_table][0][0]; for (i=0; i maxother) maxother=preprocess_like[i*(POSNUM +1) + 7]; /* This is preprocess_like[i][7] */ switch (main_method) { case STOCK_PAIR: if (mode & VER_MODE) { NEWCOILSScore(mode, sequence, &maxother,junk_scores, log_offset_to_use); /* just dump */ /* scores in junk_scores */ log_output_ver(SC2SEQ | SCSTOCK, maxscore, maxother, maxscore,seqnum, sequence.title, sequence.code, flog); } break; case MultiCoil: /****/ if (fout_coils) /*** For pos files coils, does a log file version of coil/seq score. **/ output_pos_scores(fout_coils, sequence, all_scores, log_bound, mode, by_coil_or_seq, weighted_avg); if (mode & PRN_MODE) if (Coil_Score && flog) { /*** Were having trouble with this in printing out coils, since offset *** would sometime change at ends of coils between residue offset and *** coil score offset. To make log file print out nice, use res. offset. char all_offs[MAX_NUM_SCORE_DIM][MAXSEQLEN]; get_offsets(all_scores, sequence,all_offs,3); ***********/ log_coil_output(all_scores, all_offsets,sequence, log_bound, flog, seqnum,mode); } else if (flog) log_output2_prn(mode,maxscore, maxscore, scores, scores, sequence.seqlen, seqnum, sequence, log_bound, log_bound, log_bound, flog, avg_max); /****/ break; case PairCoil: case NEWCOILS: /**** if (mode & VER_MODE) log_output_ver(SC2SEQ , maxscore, 0, maxscore,seqnum, sequence.title, sequence.code, flog); ***/ if (mode & PRN_MODE) log_output2_prn(mode,maxscore, maxscore, scores, scores,sequence.seqlen, seqnum, sequence, log_bound, log_bound, log_bound, flog, avg_max); if (flog_coil_conflicts) coil_conflicts(mode,all_scores,all_scores,sequence, seqnum,log_bound, log_bound,number_tables, flog_coil_conflicts,avg_max); break; /* case PairCoilDiff: */ case PairCoilDiffAvg: if (mode & PRN_MODE) { log_output2_prn(mode,maxscore, maxother, scores,preprocess_like, sequence.seqlen, seqnum, sequence, .5, .5, log_bound, flog, avg_max); } if (flog_coil_conflicts) { coil_conflicts(mode,all_scores,all_preprocess_like,sequence, seqnum, .5, log_bound, number_tables, flog_coil_conflicts, avg_max); } break; case ActualCoil: break; } /*************** fprintf(stderr,"\n About to go into fout."); ***************/ if ( fout) if ( (!(mode & RAW_OUT)) || ( (main_method != MultiCoil) && (main_method != PairCoilDiff))) { if (main_method== MultiCoil) tuple_output2(sequence, mode, fout, all_scores,3, 1,log_bound,3); /* Output based on 3 tables, meaning */ /* output the dimer and trimer likes */ /* and total like. */ else txt_output2(sequence, mode, fout, scores, log_bound); } else { /* Want raw tuple_output for multi_coil and PairCoilDiff*/ if (main_method == MultiCoil) if (log_offset_to_use==7) /* Combo offset just output max score */ tuple_output2(sequence, mode, fout, /** for each dimension. */ all_scores, number_tables, number_multi_lib,log_bound, number_score_dim); else /** Max reg for raw score can't really be determined. So instead */ /* for coils, output scores for the known correct register, and */ /* for pdb- output ALL possible 7 register scores (all negs). */ tuple_output_for_max(sequence,mode,fout, all_scores, number_tables, number_multi_lib, log_bound, number_score_dim); else if (main_method == PairCoilDiff) /* Output max score for each dimension. */ /** Had to do a hack since tuple_out2 expects number_lib */ /** to be array of size num_tab, so pass in 1 for num_tab */ tuple_output2(sequence, mode, fout, all_scores, 1, &number_differ_lib,log_bound, number_differ_score_dim); } } int copy_offsets(char preproc_offsets[MAXSEQLEN], Sequence sequence, char all_offsets[MAX_NUM_SCORE_DIM][MAXSEQLEN], int number_dimens, double all_scores[MAX_NUM_SCORE_DIM][MAXSEQLEN][POSNUM+1], int offset_to_use) { int tablenum,i; for (tablenum =0; tablenum < number_dimens; tablenum++) for(i=0; i max_off_score) { max_off_score=all_likes[tablenum][i][j]; max_offset=(j-i) % POSNUM; if (max_offset <0) max_offset +=POSNUM; } } if (!got_offset) offsets[tablenum][i] = max_offset; /* This is a new thing so that if an offset change is not real */ /* we keep the old offset. **********/ if ( (i>0) && (all_likes[tablenum][i][(i+offsets[tablenum][i-1]) %POSNUM] == all_likes[tablenum][i][(i+offsets[tablenum][i])%POSNUM]) ) offsets[tablenum][i] = offsets[tablenum][i-1]; } } void preprocessor_score(char lib[MAXFUNCTNUM], char multi_lib[MAX_TABLE_NUMBER][MAXFUNCTNUM], double *m, double *b, double *m_single, double *b_single, int mode, int Preprocessor_Method, int Offset_to_Use, double all_preprocess_like[MAX_NUM_SCORE_DIM][MAXSEQLEN][POSNUM+1], Sequence sequence, char all_preproc_offsets[MAX_NUM_SCORE_DIM][MAXSEQLEN], double bound) { double maxscore; char offsets[MAXSEQLEN]; int tablenum; char *local_lib; int dist; int use_pairs=1; int dimension, num_multilib; for (tablenum =0; tablenum < number_tables; tablenum++) { switch_tables(tablenum); /* Set prob. to be for that table. */ if (tablenum ==0) use_pairs =1; /* Use pairs for cctb */ else if (tablenum ==1) use_pairs = (mode & MULTI_TRIMER_PAIRS); /* Do what config file said to do for trimers. **/ switch(Preprocessor_Method) { case MultiCoil: functnum = multi_functnum[tablenum]; /* For use in scscore.c **/ num_multilib= compute_num_multilib(tablenum, number_multi_lib, mode & MULTI_TRIMER_PAIRS); for (dist=0; dist1) { local_lib = &multi_lib[tablenum][dist]; } else local_lib = multi_lib[tablenum]; dimension = compute_dimension(tablenum, number_tables, dist, number_score_dim, number_multi_lib, mode & MULTI_TRIMER_PAIRS); MultiCoilDimensionScore(mode, sequence, local_lib, &maxscore,tablenum, offsets, all_preprocess_like[dimension], Offset_to_Use,use_pairs); } break; case PairCoil: if (mode & PAIRCOIL_PAIRS) { functnum = pair_functnum; /* For use in scscore.c **/ PairCoilScore(mode, sequence, lib, m[tablenum], b[tablenum], &maxscore,tablenum, offsets, all_preprocess_like[tablenum], Offset_to_Use, 0); } else /* Use singles probabilities... for likelihood use line 0 */ SingleCoilScore(mode,sequence,m_single[tablenum], b_single[tablenum], &maxscore, tablenum, all_preprocess_like[tablenum], Offset_to_Use, mode & RAW_OUT); break; case NEWCOILS: NEWCOILSScore(mode, sequence, &maxscore, all_preprocess_like[tablenum], Offset_to_Use); break; case HMMCoil: HMMScore(sequence.seq, sequence.seqlen, hmm, all_preprocess_like[tablenum], tablenum, Offset_to_Use); break; case ActualCoil: /* if (mode & POS_MODE) */ /* Score the coils in the posfile based on */ /* their having 100% likelihood. */ ActualCoils(sequence, all_preprocess_like[tablenum], Offset_to_Use,1); break; default: fprintf(stderr,"\nBad preprocessor, using PairCoil.\n"); functnum = pair_functnum; /* For use in scscore.c **/ PairCoilScore(mode, sequence, lib, m[tablenum], b[tablenum], &maxscore,tablenum, offsets, all_preprocess_like[tablenum], Offset_to_Use, 0); break; } } if (Preprocessor_Method !=MultiCoil) get_offsets(all_preprocess_like, sequence,all_preproc_offsets, number_tables); else { /***** Get copy of res scores in tablenum and tablenum+3 **/ switch_gauss_param(Offset_to_Use, Preprocessor_Method); convert_raw_scores_to_gauss_prob_like2(all_preprocess_like,sequence.seqlen, number_tables,means_submatrix,inverse_covars, determinants, number_classes, init_class_prob, number_score_dim, mode & ONLY_COILED_CLASSES, bound,Offset_to_Use, 1); get_offsets(all_preprocess_like, sequence,all_preproc_offsets,3); if (Coil_Score) { /* Value of 1 means do average, 2 means do max. */ /*** Need to save residue scores, since the coil scores **/ /*** will replace them, but old residue scores might be needed **/ /*** again in tablenum+3. ****/ for (tablenum =0; tablenum <3; tablenum++) { average_score_over_coils3(sequence, all_preprocess_like[tablenum+3], all_preprocess_like[tablenum], all_preprocess_like[2+3],/** Total Coil Like. **/ weighted_avg, offset_to_use, all_preproc_offsets[tablenum], start_coil_at_reg_shift, bound, Coil_Score -1); } } } } double *ScoreSeq (char lib[MAXFUNCTNUM], char differ_lib[MAXFUNCTNUM], char multi_lib[MAX_TABLE_NUMBER][MAXFUNCTNUM], double *m_paircoil, double *b_paircoil, double *m_single, double *b_single, int mode, int Table, int Method, int Preprocessor_Method, int Offset_to_Use, double *maxscore, int rescore_seq, int rescore_preproc, double bound) { /** rescore_preproc == 1 means need to rescore. If it or rescore_seq are */ /** -1 then means just need to recompute MultiCoil coil and seq scores. */ /** -2 means need to change from MultiCoil coil score back to res score. */ /** If it is another non-zero */ /** value, don't need to rescore it, but do need to rescore anything that */ /** uses it as a preproc. */ double max_paircoil_score; static char all_preproc_offsets[MAX_NUM_SCORE_DIM][MAXSEQLEN]; extern char all_offsets[MAX_NUM_SCORE_DIM][MAXSEQLEN]; extern double all_preprocess_like[MAX_NUM_SCORE_DIM][MAXSEQLEN][POSNUM+1]; static double all_scores[MAX_NUM_SCORE_DIM][MAXSEQLEN][POSNUM+1]; int i,j; int tablenum; double maxsc; char *local_lib; int dist; int use_pairs=1; int num_multilib, dimension; int local_number_tables; double *preproc_res_like; int scored_already =0; /**************************************************************************/ /********* PreprocLike is used in the 2-3 stranded differentiators. ****/ if ( (rescore_preproc==1) && (!ftotal_like)) { /* A hack to save computation time. */ preprocessor_score(lib, multi_lib, m_paircoil, b_paircoil, m_single, b_single, mode, Preprocessor_Method, Offset_to_Use, all_preprocess_like, sequence, all_preproc_offsets, bound); } /*****************************************************************************/ if (Preprocessor_Method == MultiCoil) { if ( ( rescore_preproc==-1) && Coil_Score) { /* Coil_Score value of 1 means do average, 2 means do max. */ /*** Need to save residue scores, since the coil scores **/ /*** will replace them, but old residue scores might be needed **/ /*** again in tablenum+3. ****/ for (tablenum =0; tablenum <3; tablenum++) average_score_over_coils3(sequence, all_preprocess_like[tablenum+3], all_preprocess_like[tablenum], all_preprocess_like[2+3],/** Total Coil Like. **/ weighted_avg, offset_to_use, all_preproc_offsets[tablenum], start_coil_at_reg_shift, bound, Coil_Score -1); if (mode & ONLY_COILED_CLASSES) /** Convert to olig state. **/ type_class_score_convert(sequence,all_preprocess_like,bound, mode & ONLY_COILED_CLASSES,number_tables, number_classes); } else if (rescore_preproc == -2) /** Switch back to res score.*/ for (tablenum =0; tablenum <3 ; tablenum++) for (j =0; j <=POSNUM; j++) for (i=0; ibound)) all_preprocess_like[tablenum][i][j]/= all_preprocess_like[2+3][i][j]; } } /****************************************************************************/ if (Method == PairCoilDiff) local_number_tables = 1; else local_number_tables= number_tables; /*** Clean up this conditional, since really the only thing that might ***/ /*** change for paircoildiffer when change preproc is the offset to use **/ /*** (i.e. all_scores[table][i][7]. ***/ if ((rescore_seq>0) || (!is_not_differentiator(Method) && rescore_preproc) ) { scored_already=1; for (tablenum=0; tablenum1) local_lib = &multi_lib[tablenum][dist]; else local_lib = multi_lib[tablenum]; dimension = compute_dimension(tablenum, number_tables, dist, number_score_dim, number_multi_lib, mode & MULTI_TRIMER_PAIRS); MultiCoilDimensionScore(mode, sequence, local_lib, &maxsc, tablenum,all_offsets[dimension], all_scores[dimension], Offset_to_Use,use_pairs); } break; case PairCoil: if (mode & PAIRCOIL_PAIRS) { functnum = pair_functnum; /* For use in scscore.c **/ PairCoilScore(mode, sequence, lib, m_paircoil[tablenum], b_paircoil[tablenum], &maxsc, tablenum,all_offsets[tablenum], all_scores[tablenum], Offset_to_Use, (mode & RAW_OUT)); } else /*Use singles probabilities... for likelihood use line number 0 */ SingleCoilScore(mode,sequence,m_single[tablenum], b_single[tablenum], &maxsc, tablenum, all_scores[tablenum], Offset_to_Use, (mode & RAW_OUT)); break; case NEWCOILS: NEWCOILSScore(mode,sequence,&maxsc, all_scores[tablenum], Offset_to_Use); break; case HMMCoil: HMMScore(sequence.seq, sequence.seqlen, hmm, all_scores[tablenum], tablenum, Offset_to_Use); break; case PairCoilDiff: functnum = differ_functnum; /* For use in scscore.c **/ for (dimension=0; dimension1) { local_lib = &differ_lib[dimension]; } else local_lib = differ_lib; PairCoilDifferDimension(mode, sequence, local_lib, &maxsc, all_preproc_offsets[preproc_table], all_preprocess_like[preproc_table], all_scores, all_offsets, dimension, Offset_to_Use, number_differ_score_dim, bound); } /*********************** PairCoilDiffer(mode, sequence, local_lib, &maxsc, all_preproc_offsets[tablenum], all_preprocess_like[preproc_table], Offset_to_Use, tablenum, all_scores[tablenum],0, weighted_avg,start_coil_at_reg_shift); *************************/ break; /************* case PairCoilDiffAvg: functnum = differ_functnum; PairCoilDiffer(mode, sequence, differ_lib, &maxsc, all_preproc_offsets[preproc_table], all_preprocess_like[preproc_table], Offset_to_Use, tablenum, all_scores[tablenum], 1, weighted_avg,start_coil_at_reg_shift); break; ****************/ case ActualCoil: ActualCoils(sequence, all_scores[tablenum], Offset_to_Use,0); break; } if (tablenum == Table) *maxscore = maxsc; } /*********************************************************************/ /*********************************************************************/ /********* This ends the scoring along all and all dimensions ********/ /********* Now convert the raw scores that haven't been converted ****/ /********* yet into likelihoods (i.e. the multidimensional scores.****/ /***** Note: Can use same "init_class_prob" here, because since set number **/ /*** of classes to 2, the 3rd class (Pdb-) will be ignored, so just get ****/ /*** ratio of two to three stranded likelihoods, BOTH weighted by their ****/ /*** init probaabiities. ********/ if ( (Method != MultiCoil) && (Method !=PairCoilDiff)) get_offsets(all_scores, sequence,all_offsets,number_tables); else if (!(mode & RAW_OUT)) { switch_gauss_param(Offset_to_Use, Method); /***********Do conversions on PairCoilDiff scores, using preprocessor. **/ /*********** Note for differentiator there are only 2 possible classes **/ /*********** (dimer and trimer). **********/ if (Method == PairCoilDiff) { convert_raw_scores_to_gauss_prob_like2(all_scores,sequence.seqlen, number_tables, differ_means_submatrix, differ_inverse_covars,differ_determinants, 2, init_class_prob, number_differ_score_dim, 0, bound, Offset_to_Use,0); if (Preprocessor_Method == MultiCoil) preproc_res_like = &all_preprocess_like[preproc_table+3][0][0]; else preproc_res_like = &all_preprocess_like[preproc_table][0][0]; /* Only use window scores that lie completely in coil */ if ( (which_differentiator == 3) || (which_differentiator ==4)) zero_out_non_coils(all_scores, all_preprocess_like[preproc_table], sequence.seqlen, bound); if (which_differentiator == 1) zero_out_bad_windows(differ_window_length, all_scores,all_preprocess_like[preproc_table], sequence.seqlen,bound); copy_offsets(all_preproc_offsets[preproc_table],sequence, all_offsets,2, all_scores, Offset_to_Use); if (Coil_Score) /* Value of 1 means do average, 2 means do max. */ for (tablenum=0; tablenum<2; tablenum++) { average_score_over_coils3(sequence, all_scores[tablenum+3], all_scores[tablenum], preproc_res_like, /** Total Coil Like. **/ weighted_avg, offset_to_use, all_preproc_offsets[preproc_table], start_coil_at_reg_shift, bound, Coil_Score -1); } } /**********Do conversions on MultiCoil dimension scores. ******/ else if (Method == MultiCoil) { convert_raw_scores_to_gauss_prob_like2(all_scores,sequence.seqlen, number_tables, means_submatrix,inverse_covars, determinants, number_classes,init_class_prob, number_score_dim, mode & ONLY_COILED_CLASSES, bound, Offset_to_Use,1); get_offsets(all_scores, sequence,all_offsets,3); if (Coil_Score) /* Value of 1 means do average, 2 means do max. */ for (tablenum=0; tablenum<3; tablenum++) average_score_over_coils3(sequence, all_scores[tablenum+3], all_scores[tablenum], all_scores[2+3], /** Total Coil Like. **/ weighted_avg, offset_to_use, all_offsets[tablenum], start_coil_at_reg_shift, bound, Coil_Score -1); } } } /* Ends if rescore_seq */ /*****************************************************************************/ if ( (Method == MultiCoil) || (Method == PairCoilDiff)){ if (Method == MultiCoil) if ( (rescore_seq ==1 ) || /** Signals rescored MultiCoil **/ (rescore_seq==-1)) /*Signals bound changed, need to redo seq_score */ get_seq_scores(seq_scores, sequence, all_scores[0+3], all_scores[1+3],bound, weighted_avg); if (rescore_seq==-1 && !scored_already && (Method==PairCoilDiff)) { if ( (which_differentiator == 3) || (which_differentiator ==4)) zero_out_non_coils(all_scores, all_preprocess_like[preproc_table], sequence.seqlen, bound); if (which_differentiator == 1) /* Only use window scores that lie */ zero_out_bad_windows(differ_window_length, /* completely in coil */ all_scores, all_preprocess_like[preproc_table], sequence.seqlen,bound); } if ( (rescore_seq==-1) && Coil_Score && !scored_already) { if (Method == MultiCoil) { for (tablenum=0; tablenum<3; tablenum++) { average_score_over_coils3(sequence, all_scores[tablenum+3], all_scores[tablenum], all_scores[2+3], /** Total Coil Like. **/ weighted_avg, offset_to_use, all_offsets[tablenum], start_coil_at_reg_shift, bound, Coil_Score -1); if (mode & ONLY_COILED_CLASSES) /** Convert to olig. state. **/ type_class_score_convert(sequence,all_scores,bound, mode & ONLY_COILED_CLASSES, number_tables, number_classes); } } else if (Method == PairCoilDiff) { if (Preprocessor_Method == MultiCoil) preproc_res_like = &all_preprocess_like[preproc_table+3][0][0]; else preproc_res_like = &all_preprocess_like[preproc_table][0][0]; for (tablenum=0; tablenum<2; tablenum++) average_score_over_coils3(sequence, all_scores[tablenum+3], all_scores[tablenum], preproc_res_like, /** Total Coil Like. **/ weighted_avg, offset_to_use, all_offsets[tablenum], start_coil_at_reg_shift, bound, Coil_Score -1); } } else if (rescore_seq == -2) /** Switch back to reg score.*/ for (tablenum =0; tablenum <3 ; tablenum++) for (j =0; j <=POSNUM; j++) for (i=0; ibound)) all_scores[tablenum][i][j]/= all_scores[2+3][i][j]; } } /****************************************************************************/ return (double*)all_scores; } ShowSeq(char lib[MAXFUNCTNUM], char differ_lib[MAXFUNCTNUM], char multi_lib[MAX_TABLE_NUMBER][MAXFUNCTNUM], double *m, double *b, double *m_single, double *b_single, int mode, int rescore_sequence, int rescore_preproc, double bound) { double maxscore; /* fprintf(stderr,"\nIn show seq"); */ if ( (rescore_sequence) || (rescore_preproc)) all_scores = ScoreSeq(lib, differ_lib, multi_lib, m, b,m_single, b_single, mode, table, method, preprocessor_method, offset_to_use, &maxscore, rescore_sequence, rescore_preproc, bound); /* fprintf(stderr,"\n Got all_scores in showseq"); */ scores = &all_scores[table * MAXSEQLEN*(POSNUM+1)]; /* This is all_scores[main_table][0][0] */ preprocess_like = &all_preprocess_like[preproc_table][0][0]; ShowScores(scores,&sequence,NULL, bound,mode); } /*---------------------------------------------------------------------------*/ Font font; OpenScreen(int mode, double log_bound, FILE *fin, FILE *flog, FILE* flog_coil_conflicts,FILE *fout_coils, FILE *fout, double *m, double *b, double *m_single, double *b_single, char lib[MAXFUNCTNUM], char differ_lib[MAXFUNCTNUM], char multi_lib[MAX_TABLE_NUMBER][MAXFUNCTNUM], int main_method, int main_preprocessor_method, int *seqcnt, int which_priors, int which_priorp, double prior_freq_single[MAX_TABLE_NUMBER], double prior_freq_pair[MAX_TABLE_NUMBER], int good_turing_fixed_disc, int structural_pos[POSNUM+1]) { int WIDTH, HEIGHT; XGCValues values; XSetWindowAttributes attr; Pixmap icon; static unsigned char dashl[2] = {1,1}; static unsigned char dashd[2] = {2,2}; int i; WIDTH = 597; HEIGHT = 729; display = XOpenDisplay(""); /* Check mode to see if want GUI. */ if (!display || (mode & NO_GUI)) { /** COmplete logfile, etc. if can't open xdisplay. */ if (!(mode & NO_GUI)) fprintf(stderr,"Sorry, bad display.\nMaybe your DISPLAY environment variable isn't set properly,\nor you haven't xhosted this machine.\n"); if ( fout || ( (mode & PRN_MODE) || (mode & VER_MODE)) ) { fprintf(stderr,"\nThe log file will still be written by:\n\n"); for (i = 0; motd[i][0]; ++i) /* Print out title. */ fprintf(stderr," %s\n",motd[i]); seqnum=sequp; finish_log(mode,log_bound,fin,flog,flog_coil_conflicts,fout_coils, fout, m,b, m_single, b_single, lib,differ_lib,multi_lib, main_method, main_table, main_preprocessor_method, seqcnt, avg_max, which_priors, which_priorp, prior_freq_single, prior_freq_pair, good_turing_fixed_disc, structural_pos); } sclose(fin); if (flog) sclose(flog); if (flog_coil_conflicts) sclose(flog_coil_conflicts); if (fout) sclose(fout); if (ftotal_like) sclose(ftotal_like); if (fout_coils) sclose(fout_coils); exit(0); } font = XLoadFont(display,"9x15"); window = XCreateSimpleWindow(display, RootWindow(display, 0), 0, 0, WIDTH, HEIGHT/4 - 20, 0, BlackPixel(display, 0), WhitePixel(display, 0)); help = XCreateSimpleWindow(display, RootWindow(display, 0), 0, 0, WIDTH / 5 * 4 +60, HEIGHT / 5 * 3 + 200, 0, BlackPixel(display, 0), WhitePixel(display, 0)); helpb = XCreateSimpleWindow(display, window, 90, 0, 30, 15, 1, BlackPixel(display, 0), WhitePixel(display, 0)); next = XCreateSimpleWindow(display, window, 125, 0, 30, 15, 1, BlackPixel(display, 0), WhitePixel(display, 0)); prev = XCreateSimpleWindow(display, window, 160, 0, 30, 15, 1, BlackPixel(display, 0), WhitePixel(display, 0)); autob = XCreateSimpleWindow(display, window, 195, 0, 30, 15, 1, BlackPixel(display, 0), WhitePixel(display, 0)); zoomb = XCreateSimpleWindow(display, window, 230, 0, 30, 15, 1, BlackPixel(display, 0), WhitePixel(display, 0)); printb = XCreateSimpleWindow(display, window, 265, 0, 36, 15, 1, BlackPixel(display, 0), WhitePixel(display, 0)); methodb=XCreateSimpleWindow(display, window, 306, 0, 42, 15, 1, BlackPixel(display, 0), WhitePixel(display, 0)); /*** tableb=XCreateSimpleWindow(display, window, 353, 0, 36, 15, 1, BlackPixel(display, 0), WhitePixel(display, 0)); ***/ paramb=XCreateSimpleWindow(display, window, 353, 0, 36, 30, 1, BlackPixel(display, 0), WhitePixel(display, 0)); quit = XCreateSimpleWindow(display, window, 393, 0, 30, 15, 1, BlackPixel(display, 0), WhitePixel(display, 0)); /** To get size: **/ /** Printing " trimer x.xx" takes 12 char = 12*7 =84 **/ /** Each line takes 15 of y coord. giving 3*15 =45 + 10 for blank line **/ seqscoreb=XCreateSimpleWindow(display, window, XCOORD_SEQSCORES, YCOORD_SEQSCORES, 84, 55, 1, BlackPixel(display, 0), WhitePixel(display, 0)); zoom = XCreateSimpleWindow(display, RootWindow(display, 0), 0, 0, 400, 169, 0, BlackPixel(display, 0), WhitePixel(display, 0)); XSelectInput(display, helpb, ButtonPressMask); XSelectInput(display, next, ButtonPressMask); XSelectInput(display, prev, ButtonPressMask); XSelectInput(display, autob, ButtonPressMask); XSelectInput(display, zoomb, ButtonPressMask); XSelectInput(display, printb, ButtonPressMask); XSelectInput(display, methodb, ButtonPressMask); /*** XSelectInput(display, tableb, ButtonPressMask); ****/ XSelectInput(display, paramb, ButtonPressMask); XSelectInput(display, quit, ButtonPressMask); XSelectInput(display, window, ExposureMask); XSelectInput(display, zoom, ExposureMask | KeyPressMask | ButtonPressMask); XSelectInput(display, help, ExposureMask | KeyPressMask | ButtonPressMask); attr.bit_gravity = NorthWestGravity; XChangeWindowAttributes(display,window,CWBitGravity,&attr); XMapWindow(display, window); XFlush(display); values.foreground = BlackPixel(display, 0); values.background = WhitePixel(display, 0); values.font = font; gc = XCreateGC(display, window, GCForeground | GCBackground | GCFont, &values); values.line_style = LineDoubleDash; values.font = XLoadFont(display,"6x13"); gcgray = XCreateGC(display, window, GCForeground | GCBackground | GCLineStyle | GCFont, &values); XSetDashes(display, gcgray, 1, dashl, 2); values.function = GXxor; values.foreground = 1; gcdash = XCreateGC(display, window, GCForeground | GCLineStyle | GCFunction, &values); XSetDashes(display, gcdash, 0, dashd, 2); /* icon = XCreateBitmapFromData(display,window,bart_bits,bart_width,bart_height);*/ /* XSetStandardProperties(display,window,"Coiled Coil Finder","xcoil",icon,NULL,0,NULL);*/ /* XSetStandardProperties(display,zoom,"XCoil Zoomer","xcoil",icon,NULL,0,NULL);*/ /* XSetStandardProperties(display,help,"XCoil Help","xcoil",icon,NULL,0,NULL);*/ XFlush(display); } PSDefs (FILE *ps, int seqlen) { double yscaler=1; double xscaler; int number_lines; int nominallength= 300; /* Default value for number residues per line. */ if ( (seqlen%nominallength <= 50) && (seqlen >= nominallength)) { nominallength = 350; } /*************************/ if (ps_res_per_line != 0) /* 0 signifies do the default as above. */ if (ps_res_per_line == -1) /* -1 signifies put the seq all on 1 line. */ nominallength = seqlen; else nominallength = ps_res_per_line; /* How many res per line. */ /************************/ if (seqlen ==0) number_lines=1; else number_lines = (int)(seqlen-1)/nominallength +1; /* nominallength amino acids per line. */ /* C rounds down and we want 0-nominallength to go to 1 */ /* nominallength+1 to 2*noiminallength to go to 2.. */ fputs("%!\n",ps); /*****Scale and translate to fit the page if want it to rotate page. */ /* unscaled there are 9 lines per page. */ /* .89 scaling if 300 res per line. xscaler = .89 * 300/nominallength; fputs("90 rotate\n",ps); fputs("78 -60 translate\n",ps); if (number_lines > 9) yscaler = (double) 9/number_lines; fprintf(ps,"%lf %lf scale\n", xscaler, yscaler); */ /*****Scale and translate to fit the page if don't want to rotate page. */ /* unscaled there can fit 12 lines per page, so 3600. */ xscaler = .75 * 300/nominallength; /* .75 scaling if 300 res per line. */ fputs("25 740 translate\n",ps); if (number_lines > 12) yscaler = (double) 12/number_lines; fprintf(ps,"%lf %lf scale\n", xscaler,yscaler); fputs("/black {0 setgray} def\n",ps); fputs("/gray {0.6 setgray} def\n",ps); fputs("/tt /Courier findfont 15 scalefont def\n",ps); fputs("/rm /Times-Roman findfont 12 scalefont def\n",ps); fputs("/cshow {dup stringwidth pop -2 div 0 rmoveto show}def\n",ps); fputs("/rshow {dup stringwidth pop neg 0 rmoveto show}def\n",ps); fputs("/m {neg moveto} def\n",ps); fputs("/l {neg lineto stroke} def\n",ps); fputs("/rl {neg rlineto} def\n",ps); fputs("/rl_stroke {neg rlineto stroke} def\n",ps); fputs("/bar {dup 0.0019 lt {pop} {0 exch rlineto stroke} ifelse} def\n",ps); } /*** static int xmin = 60; ***** Made it global variable *****/ /*Originally was 29 */ /** xmin determines where in window start seq. **/ static int nominallength; static int length; static int ybase; static int ydelta=57; /* If click at x,y, determine which residue position it corresponds to. Return 0 if none. */ XYtoPos (int x, int y) { int i; if (x= nominallength) return 0; for (y -= ybase+ydelta; y>0; y-=ydelta, i+=nominallength); if (y<-30) return 0; if (i>=length) return 0; return i+1; } PostoXY (int i, int *x, int *y) { --i; *x = xmin + 2*(i%nominallength); *y = ybase+ydelta + ydelta*(i/nominallength); } psstr (char *str, int len, FILE *ps) { int i; for (i=0; i 0) { /* Signals above */ sprintf(buf," Dimer %5.2lf", seq_scores[0]); /* bound. */ XDrawImageString(display,*seqscoreb,gcgray,3,35,buf,strlen(buf)); sprintf(buf," Trimer %5.2lf", seq_scores[1]); XDrawImageString(display,*seqscoreb,gcgray,3,50,buf,strlen(buf)); if (ps) { fprintf(ps,"rm setfont "); fprintf(ps,"%d %d m (Seq Prob for) show\n", 15,60 +ybase-30); fprintf(ps,"%d %d m (bound %5.2lf:) show\n", 15,70 +ybase-30,bound); fprintf(ps,"%d %d m ( Dimer %5.2lf) show\n",15,80 +ybase-30, seq_scores[0]); fprintf(ps,"%d %d m ( Trimer %5.2lf) show\n",15,90 +ybase-30, seq_scores[1]); } } else { sprintf(buf," No Residue "); XDrawImageString(display,*seqscoreb,gcgray,3,35,buf,strlen(buf)); sprintf(buf," above %5.2lf ", bound); XDrawImageString(display,*seqscoreb,gcgray,3,50,buf,strlen(buf)); if (ps) { fprintf(ps,"rm setfont "); fprintf(ps,"%d %d m (Seq Prob:) show\n", 15,60 +ybase-30); fprintf(ps,"%d %d m ( No Residue) show\n", 15, 75 +ybase-30); fprintf(ps,"%d %d m ( above %5.2lf) show\n",15,85 +ybase-30, bound); } } } draw_bound_box(double bound, Window paramb) { char buf[6]; XMapWindow(display, paramb); XDrawString(display,paramb,gcgray,3,10,"bound",5); sprintf(buf,"%5.2lf", bound); XDrawImageString(display,paramb,gcgray,3,25,buf,strlen(buf)); } ShowScores (double *scores, Sequence *sequence, FILE *ps, double bound, int mode) { int i,q, x, y0, xmax; int width; int rows; int charperline; int xcoord; char MethodTitleArray[500]; /* For printing the methodname to ps file. */ char *MethodTitle = MethodTitleArray; char PreprocTitleArray[500]; /* For printing the preproc.name to ps file. */ char *PreprocTitle = PreprocTitleArray; MethodTitleArray[0] = '\0'; /* Initialization. */ PreprocTitleArray[0]='\0'; nominallength = 300; length = sequence->seqlen; /* This makes display look nice (so don't have a line with < 50 residues. */ if ((length%nominallength <= 50) && (length>= nominallength)){ nominallength = 350; /* if (ps) fputs("-50 0 translate\n",ps); */ } /*************************/ if (ps) /* Do this stuff only for postscript printout. */ if (ps_res_per_line != 0) /* 0 signifies do the default as above. */ if (ps_res_per_line == -1) /* -1 signifies put the seq all on 1 line. */ nominallength = sequence->seqlen; else nominallength = ps_res_per_line; /* How many res per line. */ /************************/ charperline = (xmin+2*nominallength-5)/9; /*Each char is "9" units witdth */ rows = (length-1)/nominallength+1; /* Number of rows to display. */ if (rows==1) xmax = 2*length+xmin; /* The max x position in current row. */ else xmax = 2*nominallength+xmin; width = xmin+2*nominallength+ WIDTH_ADD_ON; /* Width of window to display. */ /* The +WIDTH_ADD_ON used to be +11, but */ /* was cutting of some numbers. */ x = strlen(sequence->title); XClearWindow(display,window); /**********This is where window size is adjusted to fit the sequence. ***/ /********** using width and number of rows. **************************/ XResizeWindow(display,window,width, 38+15*((x+charperline-1)/charperline) +ydelta*rows); /*************************************************************************/ XMapWindow(display, helpb); XMapWindow(display, next); XMapWindow(display, prev); XMapWindow(display, autob); XMapWindow(display, zoomb); XMapWindow(display, printb); XMapWindow(display, methodb); /**** XMapWindow(display, tableb); XMapWindow(display, paramb); Now done in draw_bound_box() ******/ XMapWindow(display, quit); XSync(display,0); { XEvent report; while (XCheckTypedEvent(display,Expose,&report)); } XDrawString(display,helpb,gcgray,3,10,"help",4); XDrawString(display,next,gcgray,3,10,"next",4); XDrawString(display,prev,gcgray,3,10,"prev",4); XDrawString(display,autob,gcgray,3,10,"auto",4); XDrawString(display,zoomb,gcgray,3,10,"zoom",4); XDrawString(display,printb,gcgray,3,10,"print",5); XDrawString(display,methodb,gcgray,3,10,"method",6); /* XDrawString(display,tableb,gcgray,3,10,"table",5); */ draw_bound_box(bound, paramb); /** draw_seqscore_box(bound, seq_scores,?? &seqscoreb, ps, ybase); MOVED TO LATER (AFTER ybase gets adjusted at end of this function). ****/ XDrawString(display,quit,gcgray,3,10,"quit",4); /************This section prints out the method name..... ********************/ /**** Add 9 to xcoord for each character. *****/ XDrawString(display,window,gc, XCOORD_SCORED_BY,15,"Scored by ",10); make_methodtitle(MethodTitle, method, table, mode); XDrawString(display,window,gc,XCOORD_SCORED_BY +90,15,MethodTitle, strlen(MethodTitle)); XDrawString(display,window,gc, XCOORD_SCORED_BY,30,"Filter by ",10); make_methodtitle(PreprocTitle, preprocessor_method, preproc_table, mode); XDrawString(display,window,gc,XCOORD_SCORED_BY + 90,30, PreprocTitle,strlen(PreprocTitle)); /***************************************************************************/ XDrawString(display,window,gc, 5,30,sequence->code,strlen(sequence->code)); for (i=0,ybase=30; x > i + charperline; i += charperline) XDrawString(display,window,gc, 5,ybase+=15,&(sequence->title[i]),charperline); XDrawString(display,window,gc, 5,ybase+15,&(sequence->title[i]),x-i); if (ps) { fprintf(ps,"tt setfont 5 15 m (%s) show\n",sequence->code); /* fprintf(ps,"175 15 m (Viewgram by XCoil) show\n"); */ /* fprintf(ps,"430 15 m (Scored by %s) show\n",methodname[method]); */ fprintf(ps,"180 15 m (Scored by %s) show\n",MethodTitle); fprintf(ps,"180 26 m (Filter by %s) show\n",PreprocTitle); for (i=0,ybase=30; x > i + charperline; i += charperline) { fprintf(ps,"5 %d m (",ybase+=15); psstr(&(sequence->title[i]),charperline,ps); fprintf(ps,") show\n"); } fprintf(ps,"5 %d m (",ybase+15); psstr(&(sequence->title[i]),x-i,ps); fprintf(ps,") show\n"); fputs("rm setfont 2 setlinewidth\n",ps); } for (i=0,y0=ybase; ireg && method==HMMCoil) { if (ps) fputs("gray\n",ps); for (i=0,y0=ybase; ireg[i] != '.') { double y = y0-30*scores[8*i+sequence->reg[i]]; XDrawLine(display, window, gcgray, x, y0, x, (int)(y+0.5)); XDrawLine(display, window, gcgray, x+1, y0-1, x+1, (int)(y+0.5)); if (ps) fprintf(ps,"%d %d m %lf bar\n",x+1, y0, y0-y); } } if (ps) fputs("black\n",ps); } **********************/ /**********Now draw in the PreProcScore in gray every other line. ****/ if (ps) fputs("gray\n",ps); for (i=0,y0=ybase; i1) { sprintf(buf,"%d",which_table+1); /* XDrawString(display,window,gc, xcoord,15,buf,strlen(buf)); */ xcoord += 9*strlen(buf); MethodTitle= strcat(MethodTitle, buf); } } sprintf(buf,""); /* Clear out buf. */ if ((offset_to_use != -1) && (offset_to_use != 7)) sprintf(buf,"%c",'a' + offset_to_use); else if (Method == MultiCoil) if (offset_to_use == -1) { strcat(MethodTitle,"Max"); xcoord += 9*3; } else if (offset_to_use == 7) { strcat(MethodTitle,"Combo"); xcoord += 9*5; } /* XDrawString(display,window,gc, xcoord,15,buf,strlen(buf)); */ xcoord += 9*strlen(buf); MethodTitle= strcat(MethodTitle, buf); } /***************************************************************************/ int is_not_differentiator(int Method) { if ( (Method == MultiCoil) || (Method == PairCoil) || (Method==HMMCoil) || (Method==NEWCOILS) || (Method == ActualCoil) ) return (1); else return(0); } log_file_header(FILE *flog, FILE *flog_coil_conflicts, int mode, int argc, char *argv[], int avg_max, double bound, int preprocessor_method, int preproc_table, int log_offset_to_use, int start_coil_at_reg_shift, int weighted_avg, int main_method, int main_table) { int i; char MethodTitleArray[500]; /* For printing the methodname to ps file. */ char *MethodTitle = MethodTitleArray; MethodTitleArray[0] = '\0'; /* Initialization. */ make_methodtitle(MethodTitle, method,main_table,mode); if (flog) { fprintf(flog,"\nInput file is %s", argc>1?argv[1]:"-"); fprintf(flog, "\n\nMethod Used is %s", MethodTitle); fprintf(flog, "\nThe bound for the coil locater is %.2lf",bound); } if (flog_coil_conflicts) { fprintf(flog_coil_conflicts,"\nInput file is %s", argc>1?argv[1]:"-"); if (mode & TST_MODE0) fprintf(flog_coil_conflicts,"\nSubtracting sequence off from cctb before scoring it."); if (mode & TST_MODE1) fprintf(flog_coil_conflicts,"\nSubtracting sequence off from trimer table before scoring it."); } if (mode & WEIGHTED_PROBS) { if (flog_coil_conflicts) fprintf(flog_coil_conflicts,"\nStatistically weighted PairCoil score is "); if (flog) fprintf(flog,"\nStatistically weighted PairCoil score is "); } else { if (flog_coil_conflicts) fprintf(flog_coil_conflicts,"\nPairCoil score is "); if (flog) fprintf(flog,"\nPairCoil score is "); } if (avg_max ==0) { if (flog_coil_conflicts) fprintf(flog_coil_conflicts,"AVERAGE of residue scores."); if (flog) fprintf(flog,"AVERAGE of residue scores."); } else { if (flog_coil_conflicts) fprintf(flog_coil_conflicts,"MAX of residue scores."); if (flog) fprintf(flog,"MAX of residue scores."); } if (mode & USE_LIKE_LINE) { if (flog_coil_conflicts) fprintf(flog_coil_conflicts,"\nPaircoil used the likelihood line."); if (flog) fprintf(flog,"\nPaircoil used the likelihood line."); } else { if (flog_coil_conflicts) fprintf(flog_coil_conflicts,"\nPaircoil used direct prob. formula instead of likelihood line."); if (flog) fprintf(flog,"\nPaircoil used direct prob. formula instead of likelihood line."); } if (flog_coil_conflicts) fprintf(flog_coil_conflicts, "\nThe distances used by PairCoil are: "); if (flog) fprintf(flog, "\nThe distances used by PairCoil are: "); for (i=0; i0) && (offsets[i-1] != -1) && (offsets[i] != offsets[i-1]) && (all_preprocess_like[tab + number_tables*(int)dist] [i][(i+offsets[i-1]) %POSNUM] == all_preprocess_like[tab+number_tables*(int)dist][i][7]) ) offsets[i] = offsets[i-1]; } /**************************************************************************/ /****** Now compute coil scores for that method (table,dist) pair. ***/ get_raw_coil_score(sequence, all_preprocess_like[tab+number_tables*(int)dist], all_scores[tab+number_tables*(int)dist], mode, avg_max,bound,offsets); } } /************* Now output the tuples of coil scores. ****/ tuple_output2(sequence, mode, fout_coils, all_scores, number_tables, number_multi_lib,bound, number_score_dim); } switch_gauss_param(int Offset_to_Use, int Method) { if (Method == MultiCoil) { if (Offset_to_Use == 7) { /* Combo offset */ determinants= &both_determinants[0][0]; inverse_covars= &both_inverse_covars[0][0][0][0]; means_submatrix= &both_means_submatrix[0][0][0]; } else { /* max offset or some fixed offset */ determinants= &both_determinants[1][0]; inverse_covars= &both_inverse_covars[1][0][0][0]; means_submatrix= &both_means_submatrix[1][0][0]; } } else if (Method == PairCoilDiff) { determinants= &differ_determinants[0]; inverse_covars= &differ_inverse_covars[0][0][0]; means_submatrix= &differ_means_submatrix[0][0]; } } void usage(char *prog_command) { int i; char *prog_name; /** Strip off directories to get just the program name. **/ i=strlen(prog_command); while ((i>=0) && (prog_command[i] != '/') ) i--; prog_name= &prog_command[i+1]; fprintf(stderr,"\nUsage: %s [dimer, trimer, dimer-1, trimer-1, negatives, config_file] [multi_lib_dim1] [multi_lib_dim2] \n",prog_name); }