/* Output from p2c 1.21alpha-07.Dec.93, the Pascal-to-C translator */ /* From input file "malopt.p" */ #include /* malopt: 2d plot of optimal alignments Tom Schneider NCI/FCRDC Bldg 469. Room 144 P.O. Box B Frederick, MD 21702-1201 (301) 846-5581 (-5532 for messages) toms@ncifcrf.gov http://www-lmmb.ncifcrf.gov/~toms/ National Cancer Institute Laboratory of Mathematical Biology */ /* end of program */ /* begin module version */ #define version 1.09 /* of malopt.p 1996 March 22 origin 1996 March 21 */ /* end module version */ /* begin module describe.malopt */ /* name malopt: 2d plot of optimal alignments synopsis malopt(optalign: in, maloptp: in, xyin: out, output: out) files optalign: output of malign maloptp: parameters to control this program. Lines: 1: doclasses: integer number of classes to process 2: colorbound: highest distance of the color scale 3: xrectangle, yrectangle: size of rectangles of colors. Values of 1 and 1 give nice squares that just fit together. Values of 0.9 and 0.9 squares that don't touch each other. (This would allow more precise location of particluar pieces of data.) xyin: input to xyplo. Columns: 1: x class 2: y class 3: x H, uncertainty in bits 4: y H, uncertainty in bits 5: x rectangle width 6: y rectangle height 7: distance from x class to y class 8: hue for color computed from the distance 9: saturation for color 10: brightness for color A color scale is put on the left side of the graph, running from 0 to colorbound. The lower right triangle of the graph gives the distances according to these colors. output: messages to the user description How can we look at all the classes produced by malign? We can first normalize each alignment vector by subtracting the first value from all other values. Then we can ask what the distance is between all classes. These data are put into xyin for plotting by the xyplo program. The distances are computed and reported to xyin. The distances are also linearly transformed to give a color scale. The user can define the upper bound of the color scale, to be able to look at details for smaller distances. Distances larger than this "colorbound" are set to black. examples documentation see also malign.p, maloptp, xyplo.p, xyplop.malopt author Thomas Dana Schneider bugs technical notes See xyplo for notes on color. */ /* end module describe.malopt */ /* begin module malopt.const */ #define maxclasses 100 /* maximum number of alignment classes that can be handled */ #define maxsequences 300 /* maximum number of sequences that can be handled */ #define nwid 4 /* width of integers */ #define wid 10 /* width of numbers */ #define dec 5 /* decimal places of numbers */ /* end module malopt.const */ Static _TEXT optalign; /* file used by this program */ Static _TEXT maloptp, xyin; /* file used by this program */ Static jmp_buf _JL1; /* begin module halt */ Static Void halt() { /* stop the program. the procedure performs a goto to the end of the program. you must have a label: label 1; declared, and also the end of the program must have this label: 1: end. examples are in the module libraries. this is the only goto in the delila system. */ printf(" program halt.\n"); longjmp(_JL1, 1); } /* end module halt version = 4.15; (@ of prgmod.p 1994 November 12 */ /* begin module copyaline */ Static Void copyaline(fin, fout) _TEXT *fin, *fout; { /* copy a line from file fin to file fout */ while (!P_eoln(fin->f)) { putc(P_peek(fin->f), fout->f); getc(fin->f); } fscanf(fin->f, "%*[^\n]"); getc(fin->f); putc('\n', fout->f); } /* copyaline */ /* Local variables for themain: */ struct LOC_themain { _TEXT *xyin; long brightness; /* color brightness */ double colorbound; /* highest distance to give colors for */ double distance; /* a distance between two classes */ double H[maxclasses + 1]; /* uncertainty for each class */ long saturation; /* color saturation */ /* the vectors of each class */ long xclass; /* x class of alignment we are working on */ double xrectangle; /* the length of the x axis of the rectangles to be drawn by xyplo */ long yclass; /* y class of alignment we are working on */ double yrectangle; } ; Local double hue(color, LINK) double color; struct LOC_themain *LINK; { /* color hue. User linear transform gives range to use. Then the next transform gives a good color range. */ return (color / LINK->colorbound); } Local Void columndefinition(LINK) struct LOC_themain *LINK; { /* define xyin columns */ fprintf(LINK->xyin->f, "*\n"); fprintf(LINK->xyin->f, "* data columns:\n"); fprintf(LINK->xyin->f, "* 1: x class\n"); fprintf(LINK->xyin->f, "* 2: y class\n"); fprintf(LINK->xyin->f, "* 3: x H, uncertainty in bits\n"); fprintf(LINK->xyin->f, "* 4: y H, uncertainty in bits\n"); fprintf(LINK->xyin->f, "* 5: x rectangle width\n"); fprintf(LINK->xyin->f, "* 6: y rectangle height\n"); fprintf(LINK->xyin->f, "* 7: distance from x class to y class\n"); fprintf(LINK->xyin->f, "* 8: distance converted to color\n"); fprintf(LINK->xyin->f, "* 9: hue for color\n"); fprintf(LINK->xyin->f, "* 10: saturation for color\n"); fprintf(LINK->xyin->f, "* 11: brightness for color\n"); fprintf(LINK->xyin->f, "*\n"); } Local Void generateoutput(LINK) struct LOC_themain *LINK; { /* generate the output to xyin */ /* compute cutoff for colors */ if (LINK->distance > LINK->colorbound) LINK->brightness = 0; else LINK->brightness = 1; fprintf(LINK->xyin->f, " %*ld", nwid, LINK->xclass); fprintf(LINK->xyin->f, " %*ld", nwid, LINK->yclass); fprintf(LINK->xyin->f, " %*.*f", wid, dec, LINK->H[LINK->xclass]); fprintf(LINK->xyin->f, " %*.*f", wid, dec, LINK->H[LINK->yclass]); fprintf(LINK->xyin->f, " %*.*f", nwid, dec, LINK->xrectangle); fprintf(LINK->xyin->f, " %*.*f", nwid, dec, LINK->yrectangle); fprintf(LINK->xyin->f, " %*.*f", wid, dec, LINK->distance); fprintf(LINK->xyin->f, " %*.*f", wid, dec, hue(LINK->distance, LINK)); fprintf(LINK->xyin->f, " %*ld", nwid, LINK->saturation); fprintf(LINK->xyin->f, " %*ld\n", nwid, LINK->brightness); } /* end module copyaline version = 4.15; (@ of prgmod.p 1994 November 12 */ /* begin module malopt.themain */ Static Void themain(optalign, maloptp, xyin_) _TEXT *optalign, *maloptp, *xyin_; { /* the main procedure of the program */ struct LOC_themain V; long avalue; /* a value of the vector */ long class_; /* supposedly the class we are working on */ long classes = 0; /* number of classes of alignments */ double d; /* square of difference */ long doclasses; /* number of classes of alignments to handle */ long s; /* index to sequences */ long sequences; /* number of sequences being aligned */ double sumsq; /* sum of squares */ long theclasses[maxclasses][maxsequences]; /* the length of the y axis of the rectangles to be drawn by xyplo */ long zerovalue; /* the first value of the vector for normalizing */ long FORLIM1; V.xyin = xyin_; printf("malopt %4.2f\n", version); if (*optalign->name != '\0') { if (optalign->f != NULL) optalign->f = freopen(optalign->name, "r", optalign->f); else optalign->f = fopen(optalign->name, "r"); } else rewind(optalign->f); if (optalign->f == NULL) _EscIO2(FileNotFound, optalign->name); RESETBUF(optalign->f, Char); /* read parameters */ if (*maloptp->name != '\0') { if (maloptp->f != NULL) maloptp->f = freopen(maloptp->name, "r", maloptp->f); else maloptp->f = fopen(maloptp->name, "r"); } else rewind(maloptp->f); if (maloptp->f == NULL) _EscIO2(FileNotFound, maloptp->name); RESETBUF(maloptp->f, Char); fscanf(maloptp->f, "%ld%*[^\n]", &doclasses); getc(maloptp->f); fscanf(maloptp->f, "%lg%*[^\n]", &V.colorbound); getc(maloptp->f); fscanf(maloptp->f, "%lg%lg%*[^\n]", &V.xrectangle, &V.yrectangle); getc(maloptp->f); /* determine number of classes */ while (P_peek(optalign->f) == '*') { fscanf(optalign->f, "%*[^\n]"); getc(optalign->f); } fscanf(optalign->f, "%*[^\n]"); getc(optalign->f); /* skip blank line */ fscanf(optalign->f, "%*[^\n]"); getc(optalign->f); /* skip sequences */ fscanf(optalign->f, "%*[^\n]"); getc(optalign->f); /* skip max R line */ while (!BUFEOF(optalign->f)) { /* count blank lines */ if (P_eoln(optalign->f)) classes++; fscanf(optalign->f, "%*[^\n]"); getc(optalign->f); } /* start output */ if (*V.xyin->name != '\0') { if (V.xyin->f != NULL) V.xyin->f = freopen(V.xyin->name, "w", V.xyin->f); else V.xyin->f = fopen(V.xyin->name, "w"); } else { if (V.xyin->f != NULL) rewind(V.xyin->f); else V.xyin->f = tmpfile(); } if (V.xyin->f == NULL) _EscIO2(FileNotFound, V.xyin->name); SETUPBUF(V.xyin->f, Char); fprintf(V.xyin->f, "* malopt %4.2f\n", version); if (*optalign->name != '\0') { if (optalign->f != NULL) optalign->f = freopen(optalign->name, "r", optalign->f); else optalign->f = fopen(optalign->name, "r"); } else rewind(optalign->f); if (optalign->f == NULL) _EscIO2(FileNotFound, optalign->name); RESETBUF(optalign->f, Char); while (P_peek(optalign->f) == '*') copyaline(optalign, V.xyin); fscanf(optalign->f, "%*[^\n]"); getc(optalign->f); /* skip blank line */ fscanf(optalign->f, "%ld%*[^\n]", &sequences); getc(optalign->f); fscanf(optalign->f, "%*[^\n]"); getc(optalign->f); /* skip max R line */ fprintf(V.xyin->f, "* %ld classes found\n", classes); printf("%ld classes found\n", classes); fprintf(V.xyin->f, "* %ld sequences aligned\n", sequences); printf("%ld sequences aligned\n", sequences); fprintf(V.xyin->f, "*\n"); fprintf(V.xyin->f, "* maloptp:\n"); fprintf(V.xyin->f, "* %10ld1: doclasses: integer number of classes to process\n", doclasses); fprintf(V.xyin->f, "* %*.*f2: colorbound: maximum color distance\n", wid, dec, V.colorbound); columndefinition(&V); if (doclasses > classes) { printf("You asked to process %ld but there are only %ld classes\n", doclasses, classes); printf("%ld will be processed\n", classes); doclasses = classes; } if (doclasses > maxclasses) { printf("You asked to process %ld but I can only handle %ld classes\n", doclasses, (long)maxclasses); printf("%ld will be processed\n", (long)maxclasses); doclasses = maxclasses; } /* generate color scale */ V.saturation = 1; V.brightness = 1; for (V.yclass = 0; V.yclass < doclasses; V.yclass++) { V.xclass = 0; V.distance = V.colorbound * ((double)V.yclass / doclasses); /* set these into the H columns in case those are used */ V.H[V.xclass] = 0.0; V.H[V.yclass] = V.yclass; generateoutput(&V); } /* load up class vectors */ for (V.xclass = 1; V.xclass <= doclasses; V.xclass++) { fscanf(optalign->f, "%*[^\n]"); getc(optalign->f); /* skip blank line */ fscanf(optalign->f, "%ld", &class_); /* locate the = sign: */ while (P_peek(optalign->f) != '=') getc(optalign->f); getc(optalign->f); /* move past the '=' */ fscanf(optalign->f, "%lg%*[^\n]", &V.H[V.xclass]); getc(optalign->f); if (class_ != V.xclass) { printf("class misalignment\n"); halt(); } /* read the value for zeroing */ fscanf(optalign->f, "%ld", &zerovalue); theclasses[V.xclass-1][0] = 0; /* subtract zerovalue from itself! */ for (s = 1; s < sequences; s++) { fscanf(optalign->f, "%ld", &avalue); theclasses[V.xclass-1][s] = avalue - zerovalue; } fscanf(optalign->f, "%*[^\n]"); getc(optalign->f); } printf("class vectors loaded\n"); for (V.xclass = 1; V.xclass <= doclasses; V.xclass++) { FORLIM1 = V.xclass; for (V.yclass = 1; V.yclass <= FORLIM1; V.yclass++) { /* compute distance */ sumsq = 0.0; for (s = 0; s < sequences; s++) { d = theclasses[V.xclass-1][s] - theclasses[V.yclass-1][s]; sumsq += d * d; } V.distance = sqrt(sumsq); /* compute cutoff for colors */ if (V.distance > V.colorbound) V.brightness = 0; else V.brightness = 1; generateoutput(&V); } } printf("done\n"); } /* end module malopt.themain */ main(argc, argv) int argc; Char *argv[]; { PASCAL_MAIN(argc, argv); if (setjmp(_JL1)) goto _L1; xyin.f = NULL; strcpy(xyin.name, "xyin"); maloptp.f = NULL; strcpy(maloptp.name, "maloptp"); optalign.f = NULL; strcpy(optalign.name, "optalign"); themain(&optalign, &maloptp, &xyin); _L1: if (optalign.f != NULL) fclose(optalign.f); if (maloptp.f != NULL) fclose(maloptp.f); if (xyin.f != NULL) fclose(xyin.f); exit(EXIT_SUCCESS); } /* End. */