--- /dev/null
+/* Ethan Wolf 1995 */
+
+#include <stdio.h>
+#include <math.h>
+#include "scatter_ps.h"
+
+/*** To print a point do ". NEWCOIL_SCORE PAIRCOIL_SCORE\n" **/
+/*** Then end list of all points with "showpage"**/
+/*** The parameters are:
+ PairCoil likelihood linear region has equation like(x)=ax+b
+ where the linear region runs from likelihood lower_like to upper_like.
+
+ For computing NewCoil likelihoods we use gaussian approximations for their
+ scores on the negatives in genbnk with parameter mean-g and sd-g
+ and for the positives with parameter mean-cc and sd-cc.
+ Then their likelihood is gauss-cc(x)/[gauss-cc(x) + 30*gauss-g(x)]
+ (assuming that negatives are 30 times more prevalent than coils.
+******/
+
+void scatter_trailer(FILE *scatter)
+{
+ fprintf(scatter,"showpage\n");
+}
+
+
+
+
+/** Read in scores output by "scboth -v" which are in format:
+ "seq# stock_score paircoil_score seq_code seq_name".
+ Then output them as points in scatter plot as "stock_score paircoil_score ."
+**/
+void scatter_points(FILE *scatter, FILE *scores)
+{
+ char line[300];
+ double stock_score, paircoil_score;
+
+ while (fgets(line,300,scores)) {
+ sscanf(line,"%*d %lf %lf", &stock_score, &paircoil_score);
+ fprintf(scatter,"%lf %lf .\n", stock_score, paircoil_score);
+ }
+
+}
+
+
+void scatter_header(FILE *scatter, double a, double b, double lower_like, double upper_like)
+{
+ fprintf(scatter,"%%!PS-Adobe-\n");
+ fprintf(scatter,"%%%%Center and scale the figure\n");
+ fprintf(scatter,"8.5 36 mul 11 36 mul translate\n");
+ fprintf(scatter,"0.8 0.8 scale\n");
+ fprintf(scatter,"-8.5 36 mul -11 36 mul translate\n");
+
+ fprintf(scatter,"\n%%%%The Figure Title.\n");
+ fprintf(scatter,"/name () def\n");
+ fprintf(scatter,"\n%%For printing on vertical and horizontal axis.\n");
+ fprintf(scatter,"/rshow {dup stringwidth pop neg 0 rmoveto show} def\n");
+ fprintf(scatter,"/cshow {dup stringwidth pop -2 div 0 rmoveto show} def\n");
+ fprintf(scatter,"/label {/Times-Roman findfont 20 scalefont setfont cshow} def\n");
+ fprintf(scatter,"/ttlabel {/Courier findfont 20 scalefont setfont cshow} def\n");
+ fprintf(scatter,"%%%%Center and print the Title.\n");
+ fprintf(scatter,"4.5 72 mul 9.5 72 mul moveto name label\n");
+
+ fprintf(scatter,"\n/adjust {ymin sub 10 mul exch xmin sub 10 mul exch} def\n");
+
+
+ fprintf(scatter,"\n%%%%The following is used for printing out the scatter plot points.\n");
+/** Unnecessary Stuff??
+ fprintf(scatter,"/choose { pop true } def %% choose will be redefined for selective plotting\n");
+ fprintf(scatter,"/. {choose {adjust .25 0 360 arc fill} {pop pop true} ifelse} def\n");
+**/
+ fprintf(scatter,"/. {1 index 1 index ymin gt exch xmin gt and {adjust .9 0 360 arc fill} {pop pop}ifelse} def\n\n");
+
+
+ fprintf(scatter,"/xmin -2 def\n");
+ fprintf(scatter,"/ymin -26 def\n");
+ fprintf(scatter,"/xmax 30 def\n");
+ fprintf(scatter,"/ymax 28 def\n");
+ fprintf(scatter,"/vline 1.35 ln 28 mul def\n");
+
+ fprintf(scatter,"/e 2.7182818284590452353 def\n\n");
+
+ fprintf(scatter,"%%%%For computing down*gaussian(x) with input x mean sd\n");
+ fprintf(scatter,"/tgauss { \n");
+ fprintf(scatter," /stddev exch def sub stddev div dup mul -2 div %% -[(x-mean)/sd]^2/2\n");
+ fprintf(scatter," e exch exp %% e^<above>\n");
+ fprintf(scatter," 2.5066282746310002 div stddev div down mul\n");
+ fprintf(scatter,"} def\n\n");
+
+
+ stock_likelihood_setup(scatter);
+
+
+ define_paircoil_likelihoodline(scatter, a, b, lower_like, upper_like);
+
+
+
+ fprintf(scatter,"\n/axes { %% Do the axes\n");
+ fprintf(scatter," 90 100 translate\n");
+ fprintf(scatter,"%%%% Horizontal\n");
+ fprintf(scatter," 0.5 setlinewidth\n");
+ fprintf(scatter," xmin ymin adjust moveto xmax ymin adjust lineto stroke\n");
+ fprintf(scatter," xmax ymin adjust moveto -8 4 rlineto stroke\n");
+ fprintf(scatter," xmax ymin adjust moveto -8 -4 rlineto stroke\n");
+ fprintf(scatter," vline ymin adjust moveto 0 -66 rmoveto\n");
+ fprintf(scatter," (NEWCOILS) ttlabel\n");
+ fprintf(scatter,"%%%% Vertical\n");
+ fprintf(scatter," xmin ymin adjust moveto xmin ymax adjust lineto stroke\n");
+ fprintf(scatter," xmin ymax adjust moveto 4 -8 rlineto stroke\n");
+ fprintf(scatter," xmin ymax adjust moveto -4 -8 rlineto stroke\n");
+ fprintf(scatter," xmin hline adjust moveto -80 0 rmoveto\n");
+ fprintf(scatter," gsave 90 rotate (PairCoil) ttlabel grestore\n");
+ fprintf(scatter,"%%%% Vertical Ticks\n");
+ fprintf(scatter," -24 2 26 {xmin exch adjust moveto -3 0 rmoveto 6 0 rlineto stroke} for\n");
+ fprintf(scatter,"%%%% Horizontal Ticks\n");
+ fprintf(scatter," [10 11 12 13 14 15 16 18 20 22 24 26 28]\n");
+ fprintf(scatter," {10 div ln 28 mul ymin adjust moveto 0 -3 rmoveto 0 6 rlineto stroke} forall\n");
+ fprintf(scatter,"%%%% Horizontal Numbers\n");
+ fprintf(scatter," /Times-Roman findfont 15 scalefont setfont\n");
+ fprintf(scatter," [10 11 12 13 14 16 18 20 24 28]\n");
+ fprintf(scatter," {10 div dup ln 28 mul ymin adjust -22 add moveto\n");
+ fprintf(scatter," 5 string cvs cshow} forall\n");
+
+
+ stock_percentages(scatter);
+
+
+ fprintf(scatter,"%%%% Vertical Numbers (paircoil scores)\n");
+ fprintf(scatter," -24 2 26 {dup xmin exch adjust moveto -10 -4 rmoveto\n");
+ fprintf(scatter," 4 string cvs rshow} for\n");
+ fprintf(scatter,"%%%% Vertical Likelihoods (percentages)\n");
+ fprintf(scatter," -10 2 18 {dup xmin exch adjust moveto -10 -4 rmoveto -60 0 rmoveto sfit show (%%) show} for\n");
+ fprintf(scatter," newpath\n");
+ fprintf(scatter,"} def %%Ends the Definition of the Axes\n\n\n");
+
+
+
+ fprintf(scatter,"0 setgray\n");
+ fprintf(scatter,"axes\n");
+ fprintf(scatter,"newpath\n");
+ fprintf(scatter,"vline dup\n");
+ fprintf(scatter,"ymin adjust moveto ymax adjust lineto\n");
+ fprintf(scatter,"stroke\n");
+ fprintf(scatter,"hline dup\n");
+ fprintf(scatter,"xmin exch adjust moveto xmax exch adjust lineto\n");
+ fprintf(scatter,"stroke\n");
+ fprintf(scatter,"newpath\n");
+
+}
+
+
+
+
+
+/***This prints out post script stuff for the least square line ***/
+/***between 20% and 80% determined by the rest of the likelihood program ***/
+
+void define_paircoil_likelihoodline(FILE *scatter, double a, double b, double lower_like, double upper_like)
+{
+ /** x sfit should ouput the likelihood score for x in the linear **/
+ /** region of the likelihood graph (between lower_like and upper_like) **/
+ /** where the curve is like=ax+b. **/
+
+ fprintf(scatter,"\n%%%%x sfit should ouput the likelihood score for x in the linear region of the\n");
+ fprintf(scatter,"%%%%likelihood graph (between lower_like and upper_like) where the curve is like=ax+b.\n");
+
+ fprintf(scatter,"/a %lf def\n", a);
+ fprintf(scatter,"/b %lf def\n", b);
+ fprintf(scatter,"%%%%hline is used to draw the 50 percent likelihood line for PairCoil\n");
+ fprintf(scatter,"/hline .5 b sub a div def\n");
+ fprintf(scatter,"/lower_like %.0lf def\n",lower_like);
+ fprintf(scatter,"/upper_like %.0lf def\n",upper_like);
+ fprintf(scatter,"/sfit {\n");
+ fprintf(scatter," a mul b add 100 mul\n");
+ fprintf(scatter," .5 add cvi %%Rounding to nearest percentage.\n");
+ fprintf(scatter," dup upper_like gt %%If likelihood > upper_like\n");
+ fprintf(scatter," {(>) show upper_like 3 string cvs}\n");
+ fprintf(scatter," %%then print > and make string upper_like.\n");
+ fprintf(scatter," {dup lower_like lt %%else if < lower_like \n");
+ fprintf(scatter," {(<) show lower_like 3 string cvs}\n");
+ fprintf(scatter," %%then print < and make string lower_like.\n");
+ fprintf(scatter," {( ) show 3 string cvs} ifelse %%else print a space and\n");
+ fprintf(scatter," %%make likelihood a string and convert to text\n");
+ fprintf(scatter," } ifelse %%ends else for when upper_like>likelihood(x)\n");
+ fprintf(scatter,"} def\n\n");
+
+}
+
+
+void stock_likelihood_setup(FILE *scatter)
+{
+ fprintf(scatter,"%%%%To compute the likelihood of STOCK use:\n");
+ fprintf(scatter,"%%%% gauss_positives(x)/[gauss_pos(x)+30*gauss_genbnk(x)]\n");
+
+ fprintf(scatter,"%%%%NewCoil parameters for computing likelihoods.\n");
+ fprintf(scatter,"/mean-g 0.77 def\n");
+ fprintf(scatter,"/sd-g 0.20 def\n");
+ fprintf(scatter,"/mean-cc 1.63 def\n");
+ fprintf(scatter,"/sd-cc 0.24 def\n");
+ fprintf(scatter,"/down 1 def\n\n");
+
+ fprintf(scatter,"/prob { dup mean-g sd-g tgauss 30 mul exch\n");
+ fprintf(scatter," mean-cc sd-cc tgauss dup 3 1 roll add div\n");
+ fprintf(scatter,"} def\n");
+ fprintf(scatter,"/prob-log {28 div 2.7182818284590452353 exch exp prob} def\n\n");
+}
+
+
+void stock_percentages(FILE *scatter)
+{ fprintf(scatter,"%%%% Stock Percentages at 1.2, 1.3, 1.4.\n");
+ fprintf(scatter," 1.1 ln 28 mul ymin adjust -44 add moveto (< 1%) cshow\n");
+ fprintf(scatter," 1.6 ln 28 mul ymin adjust -44 add moveto (>99%) cshow\n");
+ fprintf(scatter," [12 13 14]\n");
+ fprintf(scatter," {10 div ln 28 mul dup ymin adjust -44 add moveto\n");
+ fprintf(scatter," prob-log 100 mul cvi 12 string cvs cshow} forall\n\n");
+}
git://source.jalview.org / jpred.git / blobdiff