src/hg/instinct/hgBamBam/hgAnnotations.c 1.2
1.2 2010/05/26 05:20:33 jsanborn
updated
Index: src/hg/instinct/hgBamBam/hgAnnotations.c
===================================================================
RCS file: /projects/compbio/cvsroot/kent/src/hg/instinct/hgBamBam/hgAnnotations.c,v
retrieving revision 1.1
retrieving revision 1.2
diff -b -B -U 4 -r1.1 -r1.2
--- src/hg/instinct/hgBamBam/hgAnnotations.c 25 May 2010 20:22:44 -0000 1.1
+++ src/hg/instinct/hgBamBam/hgAnnotations.c 26 May 2010 05:20:33 -0000 1.2
@@ -9,8 +9,9 @@
#define EXPR_DATA_SHADES 16
#define DEFAULT_MAX_DEVIATION 0.7
#define COLOR_SCALE 1
+#define GENE_HEIGHT 10
#include <limits.h>
#include "common.h"
#include "cart.h"
@@ -19,8 +20,9 @@
#include "hCommon.h"
#include "hdb.h"
#include "cytoBand.h"
#include "hCytoBand.h"
+#include "spaceSaver.h"
#include "hgBamBam.h"
static char const rcsid[] = "$Id$";
@@ -46,8 +48,10 @@
while ((row = sqlNextRow(sr)) != NULL)
{
struct genePred *gp = genePredLoad(row+1);
tuple = bedFromGenePred(gp);
+ if (row[12])
+ tuple->name = cloneString(row[12]); // switch to name2;
genePredFree(&gp);
slAddHead(&tupleList, tuple);
}
@@ -68,8 +72,17 @@
else
return 0;
}
+void drawBox(struct vGfx *vg, double x1, double x2, int pStart, int y, int h, Color valCol)
+{
+int x = round(x1) + pStart;
+int w = round(x2) - round(x1);
+if (w == 0)
+ w = 1;
+vgBox(vg, x, y, w, h, valCol);
+}
+
void addHmPixel(struct hmPixel **hmList, struct hash *pixelHash,
char mod, double x1, double x2, int pStart, int y, int h)
{
/* calculate starting position and width of box to draw */
@@ -259,19 +272,166 @@
vgBox(vg, hm->x, hm->y, hm->w, hm->h, valCol);
}
}
-
-char *annotationGif(struct settings *settings, char *tableName)
-/* Create genome GIF file and HT that includes it. */
+struct spaceSaver *bedToSS(struct settings *settings, char *db,
+ char *tableName, boolean drawLabels)
{
-if (!settings || !tableName)
+if (!settings)
return NULL;
-struct hvGfx *vg;
+/* get remote database connection */
+struct sqlConnection *conn = hAllocConn(db);
-char *db = "hg18";
+struct bed *nb, *ghBed = NULL;
+
+struct chromLay *cl = settings->layoutList;
+
+ghBed = getTable(conn, tableName, cl->chrom);
+if (!ghBed)
+ return NULL;
+
+struct bed *filterList =
+ bedFilterListInRange(ghBed, NULL, cl->chrom, cl->baseStart, cl->baseEnd);
+
+struct spaceSaver *ss = spaceSaverNew(cl->pxStart, cl->pxStart + cl->pxWidth, 100);
+
+double pixPerBase = (double) settings->width / (double) settings->totalBases;
+
+int baseStart = cl->baseStart;
+
+for(nb = filterList; nb; nb = nb->next)
+ {
+ int cStart = nb->chromStart;
+ int cEnd = nb->chromEnd;
+
+ int x1 = round(pixPerBase * (double) (cStart - baseStart));
+ int x2 = round(pixPerBase * (double) (cEnd - baseStart));
+
+ if (drawLabels)
+ x1 -= mgFontStringWidth(settings->font, nb->name) + 4;
+ spaceSaverAdd(ss, x1, x2, nb);
+ }
+spaceSaverFinish(ss);
+
+return ss;
+}
+
+void drawBedTable(struct vGfx *vg, struct settings *settings,
+ struct spaceSaver *ss, boolean drawLabels)
+/* Draw chromosome graph on all chromosomes in layout at given
+ * y offset and height. */
+{
+double pixPerBase = (double) settings->width / (double) settings->totalBases;
+
+Color upShades[EXPR_DATA_SHADES];
+static struct rgbColor zeroColor = {0, 0, 0}; // black
+static struct rgbColor highColor = {255, 0, 0}; // red
+vgMakeColorGradient(vg, &zeroColor, &highColor, EXPR_DATA_SHADES, upShades);
+
+struct chromLay *cl = settings->layoutList;
+
+int pStart = cl->pxStart;
+int baseStart = cl->baseStart;
+
+struct spaceNode *sn;
+for (sn = ss->nodeList; sn; sn = sn->next)
+ {
+ int row = sn->row * GENE_HEIGHT;
+ struct bed *nb = (struct bed *) sn->val;
+
+ int cStart = nb->chromStart;
+ int cEnd = nb->chromEnd;
+
+ double x1, x2;
+ int baseline = cStart - baseStart;
+
+ x1 = pixPerBase * (double) (cStart - baseStart);
+ x2 = pixPerBase * (double) (cEnd - baseStart);
+
+ if (drawLabels)
+ vgTextRight(vg, 0, row, x1-2, settings->fontHeight,
+ MG_BLACK, settings->font, nb->name);
+
+ if ( (x2 - x1) > 2.0)
+ { /* width > two pixels, attempt to draw exon/utr structure
+ * similar to genome browser */
+ drawBox(vg, x1, x2, pStart, 4 + row, 1, MG_BLACK);
+
+ int tStart = nb->thickStart - nb->chromStart;
+ int tEnd = nb->thickEnd - nb->chromStart;
+
+ int i, bStart, bEnd;
+ for (i = 0; i < nb->blockCount; i++)
+ {
+ bStart = nb->chromStarts[i];
+ bEnd = bStart + nb->blockSizes[i];
+
+ if (bStart < tStart)
+ { /* 5' utr region */
+ x1 = pixPerBase * (double) (bStart + baseline);
+ if (bEnd < tStart)
+ x2 = pixPerBase * (double) (bEnd + baseline);
+ else
+ x2 = pixPerBase * (double) (tStart + baseline);
+
+ drawBox(vg, x1, x2, pStart, 2 + row, 5, MG_BLACK);
+
+ if (bEnd > tEnd)
+ { /* coding region inside single block */
+ x1 = pixPerBase * (double) (tStart + baseline);
+ x2 = pixPerBase * (double) (tEnd + baseline);
+ drawBox(vg, x1, x2, pStart, 0 + row, 9, MG_BLACK);
+
+ x1 = pixPerBase * (double) (tEnd + baseline);
+ x2 = pixPerBase * (double) (bEnd + baseline);
+ drawBox(vg, x1, x2, pStart, 0 + row, 9, MG_BLACK);
+ }
+ else if (bEnd > tStart)
+ { /* coding region downstream of 5' utr */
+ x1 = pixPerBase * (double) (tStart + baseline);
+ x2 = pixPerBase * (double) (bEnd + baseline);
+ drawBox(vg, x1, x2, pStart, 0 + row, 9, MG_BLACK);
+ }
+ }
+ else if (bEnd > tEnd)
+ { /* 3' utr region */
+ if (bStart > tEnd)
+ x1 = pixPerBase * (double) (bStart + baseline);
+ else
+ x1 = pixPerBase * (double) (tEnd + baseline);
+ x2 = pixPerBase * (double) (bEnd + baseline);
+
+ drawBox(vg, x1, x2, pStart, 2 + row, 5, MG_BLACK);
+ if (bStart < tEnd)
+ { /* coding region upstream of 3' utr */
+ x1 = pixPerBase * (double) (bStart + baseline);
+ x2 = pixPerBase * (double) (tEnd + baseline);
+ drawBox(vg, x1, x2, pStart, 0 + row, 9, MG_BLACK);
+ }
+ }
+ else
+ { /* block is all coding region, no utr */
+ x1 = pixPerBase * (double) (bStart + baseline);
+ x2 = pixPerBase * (double) (bEnd + baseline);
+ drawBox(vg, x1, x2, pStart, 0 + row, 9, MG_BLACK);
+ }
+ }
+ }
+ else
+ { /* width < 2 pixels -- gene too dense, so just draw one large block */
+ x1 = pixPerBase * (double) (cStart - baseStart);
+ x2 = pixPerBase * (double) (cEnd - baseStart);
+ drawBox(vg, x1, x2, pStart, 0 + row, 9, MG_BLACK);
+ }
+ }
+
+spaceSaverFree(&ss);
+}
+
+char *annotationWholeGenome(struct settings *settings, char *db, char *tableName)
+{
int width = settings->width;
int height = settings->height; // was set to 10
if (width * height == 0)
@@ -283,10 +443,12 @@
off_t size = fileSize(md5Tn.forCgi);
if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG)
{
- vg = hvGfxOpenGif(width, height, md5Tn.forCgi, FALSE);
+ struct hvGfx *vg = hvGfxOpenGif(width, height, md5Tn.forCgi, FALSE);
+
drawTable(vg->vg, settings, db, tableName, height);
+
hvGfxClose(&vg);
}
char *filename = replaceChars(md5Tn.forHtml, "..", "");
@@ -292,8 +454,62 @@
char *filename = replaceChars(md5Tn.forHtml, "..", "");
return filename;
}
+char *annotationSingleGenome(struct settings *settings, char *db, char *tableName)
+{
+struct chromLay *cl = settings->layoutList;
+
+boolean drawLabels = TRUE;
+if (cl->baseEnd - cl->baseStart > 50000000)
+ drawLabels = FALSE;
+
+struct spaceSaver *ss = bedToSS(settings, db, tableName, drawLabels);
+
+int maxRow = -1;
+struct spaceNode *sn;
+for (sn = ss->nodeList; sn; sn = sn->next)
+ {
+ if (sn->row > maxRow)
+ maxRow = sn->row;
+ }
+if (maxRow == -1)
+ return NULL;
+
+settings->height = (maxRow + 1) * GENE_HEIGHT;
+
+struct tempName md5Tn;
+char *strToHash = cartSettingsString(bbPrefix, "annotationGif");
+trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
+
+off_t size = fileSize(md5Tn.forCgi);
+if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG)
+ {
+ struct hvGfx *vg = hvGfxOpenGif(settings->width, settings->height, md5Tn.forCgi, FALSE);
+
+ drawBedTable(vg->vg, settings, ss, drawLabels);
+
+ hvGfxClose(&vg);
+ }
+
+char *filename = replaceChars(md5Tn.forHtml, "..", "");
+return filename;
+}
+
+char *annotationGif(struct settings *settings, char *tableName)
+/* Create genome GIF file and HT that includes it. */
+{
+if (!settings || !tableName)
+ return NULL;
+
+char *db = "hg18";
+
+if (slCount(settings->layoutList) > 1)
+ return annotationWholeGenome(settings, db, tableName);
+else
+ return annotationSingleGenome(settings, db, tableName);
+}
+
double basesInScale(double numBases)
{
double exp = floor(log(numBases)/log(10.0));
if (exp < 0.0)