src/hg/instinct/hgHeatmap2/hgCircleMaps.c 1.7

1.7 2009/02/27 21:26:46 jsanborn
updated hgcirclemaps
Index: src/hg/instinct/hgHeatmap2/hgCircleMaps.c
===================================================================
RCS file: /projects/compbio/cvsroot/kent/src/hg/instinct/hgHeatmap2/hgCircleMaps.c,v
retrieving revision 1.6
retrieving revision 1.7
diff -b -B -U 1000000 -r1.6 -r1.7
--- src/hg/instinct/hgHeatmap2/hgCircleMaps.c	28 Jan 2009 22:44:18 -0000	1.6
+++ src/hg/instinct/hgHeatmap2/hgCircleMaps.c	27 Feb 2009 21:26:46 -0000	1.7
@@ -1,965 +1,1113 @@
 /* hgCircleMaps.c - circular heatmaps to visualize genomic high-throughput data for 
  * multiple samples in a concise image suitable for placement in a pathway diagram
  * This file is copyright 2009 J. Zachary Sanborn.  All rights reserved. */  
 
 #define EXPR_DATA_SHADES 16
 #define COLOR_SCALE 1
 #define DEFAULT_MAX_DEVIATION 0.7
 #define PI 3.14159265358979323846
 
 #include "common.h"
 #include "bed.h"
 #include "cart.h"
 #include "errCatch.h"
 #include "hash.h"
 #include "hdb.h"
 #include "hgColors.h"
 #include "hPrint.h"
 #include "vGfx.h"
 #include "hvGfx.h"
 #include "hgHeatmapLib.h"
 #include "featuresLib.h"
 #include "hgHeatmap2.h"
 #include "sortFeatures.h"
 #include "gfxPoly.h"
 #include "json.h"
 
 static char const rcsid[] = "$Id$";
 
 struct ring {
     struct ring *next;       /* next ring */
 
     char *name;              /* name of ring (optional) */
     char *dataset;           /* name of dataset */
     char *type;              /* type of data (feature, SNP, CNV, expression, etc.) */
 
     int numElements;         /* number of elements in ring */
     double minR;             /* min radius for current ring */
     double maxR;             /* max radius for current ring */
 
     double gain;             /* gain setting for ring */
     double maxDev;           /* maximum deviation for ring */
 
     double median;           /* median value for ring */
     boolean drawSummary;     /* draw summary */
 
+    boolean visible;         /* invisible rings exist, can be sorted, not displayed */
     boolean sort;            /* Sort this ring or not */
     int sortDir;             /* +1 = ASC, -1 = DESC */
     double maxVal;           /* maximum value, used for sorting */
 
     struct rgbColor *zeroColor; /* colors for ring */
     struct rgbColor *highColor; 
     struct rgbColor *lowColor;
 
     struct hash *data;       /* data for ring, keyed by sample */
 };
 
 struct circleMap {
     char *name;              /* name of circle (e.g. gene) */
     struct slName *samples;  /* global list of samples, order of samples defined 
 			      * position in ring */
     struct ring *rings;      /* data defining all rings for circle */
 
     double minR;             /* min radius for all rings */
     double maxR;             /* max radius for all rings */
     int width;               /* pix width */
 
     boolean drawSummary;     /* draw summary (pie chart) */
+    boolean subgrouped;      /* if successful subgrouping applied */
+
     struct hash *data;       /* data for center (pie chart) */
 };
 
 char *heatMapDbProfile = "localDb";
 
 double maxDev(struct genoHeatmap *gh)
 {
 if (!gh)
     return DEFAULT_MAX_DEVIATION;
 
 if (gh->expScale > 0)
     return gh->expScale;
 
 return DEFAULT_MAX_DEVIATION;
 }
 
 void vgRadiusBox(struct vGfx *vg, int x0, int y0, double rStart, double rStop, 
 		 double tStart, double tStop, Color col)
 {
+if ((rStart == rStop) || (tStart == tStop)) // zero height/thickness
+    return;
+
 struct gfxPoly *gp = gfxPolyNew();
 
 int x, y;
 double r, t, tStep;
 
 /* rStart -> rStop @ tStart -- one point at each edge */
 for (r = rStart; r <= rStop; r += (rStop - rStart))
     {
     x = floor(r * sin(tStart));
     y = floor(r * cos(tStart));
     gfxPolyAddPoint(gp, x + x0, y + y0);
     }
 
 /* tStart -> tStop @ rStop */
 tStep = 0.5 / (2.0*PI);
 for (t = tStart; t <= tStop; t += tStep)
     {
     x = floor(rStop * sin(t));
     y = floor(rStop * cos(t));
     gfxPolyAddPoint(gp, x + x0, y + y0);
     }
 
 /* rStop -> rStart @ tStop -- one point at each edge */
 for (r = rStop; r >= rStart; r -= (rStop - rStart))
     {
     x = floor(r * sin(tStop));
     y = floor(r * cos(tStop));
     gfxPolyAddPoint(gp, x + x0, y + y0);
     }
 
 /* tStop -> tStop @ rStart */
 for (t = tStop; t <= tStart; t -= tStep)
     {
     x = floor(rStart * sin(t));
     y = floor(rStart * cos(t));
     gfxPolyAddPoint(gp, x + x0, y + y0);
     }
 
 vgDrawPoly(vg, gp, col, TRUE);
 gfxPolyFree(&gp);
 }
 
 struct circleMap *newCircleMap(char *name)
 { /* initialize a new circle layout */
 struct circleMap *cm;
 AllocVar(cm);
 cm->name = cloneString(name);
 cm->samples = NULL;
 cm->rings = NULL;
 
 cm->minR = 0.0;
 cm->maxR = 0.0;
 cm->width = 0;
 cm->data = NULL;
 
 cm->drawSummary = FALSE;
+cm->subgrouped = FALSE;
 
 return cm;
 }
 
-struct ring *addRingToCircleMap(struct circleMap *cm, char *name, struct genoHeatmap *gh)
+struct ring *cloneRing(struct ring *rg)
+{
+if (!rg)
+    return NULL;
+
+struct ring *rg2;
+AllocVar(rg2);
+rg2->name = cloneString(rg->name);
+rg2->dataset = cloneString(rg->dataset);
+rg2->type = cloneString(rg->type);
+rg2->maxDev = rg->maxDev;
+rg2->gain = rg->gain;
+rg2->numElements = rg->numElements;
+rg2->minR = rg->minR;
+rg2->maxR = rg->maxR;
+rg2->data = NULL;
+
+rg2->drawSummary = rg->drawSummary;
+rg2->median = 0.0;
+
+rg2->visible = rg->visible;
+
+rg2->sort = rg->sort;   
+rg2->sortDir = rg->sortDir;
+rg2->maxVal = rg->maxVal;
+
+return rg2;
+}
+
+struct ring *addRingToCircleMap(struct circleMap *cm, char *name, struct genoHeatmap *gh, 
+				boolean isVisible)
 {  /* Initialize a ring, add it circleLay, and return reference to it */ 
 struct ring *rg;
 AllocVar(rg);
 rg->name = cloneString(name);
 rg->dataset = cloneString(gh->name);
 rg->type = cloneString(gh->platform);
 rg->maxDev = maxDev(gh);
 rg->gain = gh->gainSet;
 rg->numElements = 0;
 rg->minR = 0.0;
 rg->maxR = 0.0;
 rg->data = NULL;
 
 rg->drawSummary = FALSE;
 rg->median = 0.0;
 
+rg->visible = isVisible;
+
 rg->sort = FALSE;   
 rg->sortDir = 0;
 rg->maxVal = 0.0;
 
 slAddHead(&cm->rings, rg);
 return rg;
 }
 
 void addGeneDataToRing(struct ring *rg, struct hash *geneHash, char *gene, 
 		       struct slName *sampleList)
 {
 if (!rg || !geneHash)
     return;
 
 struct hashEl *el = hashLookup(geneHash, gene);
 if (!el)
     return;
 rg->numElements = slCount(el);  // number of probes matching 'gene'
 
 if (!rg->data)
     rg->data = hashNew(0);
 
 struct bed *nb, *nbList = NULL;
 for ( ; el; el = el->next)      // loop through probes
     {
     struct bed *nb = el->val;
     if (!nb)
 	continue;
     slAddHead(&nbList, nb);
     }
 
 int i = 0;
 struct slName *sl;
 struct slDouble *ad, *allData = NULL;
 for (sl = sampleList; sl; sl = sl->next)
     {
     struct slDouble *sd, *sdList = NULL;
     for (nb = nbList; nb; nb = nb->next)
 	{
 	double val = nb->expScores[i];
 
 	if (val > rg->maxVal)
 	    rg->maxVal = val;
 
 	/* Store copy of all data in list to determine median value */
 	ad = slDoubleNew(val);
 	slAddHead(&allData, ad);
 
 	sd = slDoubleNew(val);
 	slAddHead(&sdList, sd);
 	}
 
     hashAdd(rg->data, sl->name, sdList);
     i++;  // increment expId pointer.
     }
 
 rg->median = slDoubleMedian(allData);
 slFreeList(&allData);
 }
 
 void addFeatureDataToRing(struct sqlConnection *conn, struct ring *rg, 
 			  struct column *col, struct genoHeatmap *gh)
 {
 if (!rg)
     return;
 
 char *labTable = gh->patTable;
 char *value = gh->sampleField;
 char *key = gh->patField;
 char *db = gh->patDb; 
 
 if ((labTable == NULL) || (key == NULL) || (value == NULL) || (db==NULL))
     return; 
 
 if (!rg->data)
     rg->data = hashNew(0);
 rg->numElements = 1;
 
 char *colorScheme = col->cellColorSchemeVal(col);
 rg->type = cloneString(colorScheme);
 
 double absVal,  minVal, maxVal;
 char *minCutVal, *maxCutVal;
 double offset;
 int reverseColor;
 
 offset = atof(col->cellOffsetVal(col));
 reverseColor = atoi(col->cellColorReverseVal(col));
 
 char *minValStr = col->cellMinVal(col, conn);
 if (!minValStr)
     minVal = 0.0;
 else
     minVal = atof(minValStr) + offset;
 
 char *maxValStr = col->cellMaxVal(col, conn);
 if (!maxValStr)
     maxVal = 0.0;
 else
     maxVal = atof(maxValStr) + offset;
 
 if (reverseColor == -1)
     {
     maxVal = reverseColor*maxVal;
     minVal = reverseColor*maxVal;
     }
 minCutVal = col->cellMinCutVal(col, conn);
 maxCutVal = col->cellMaxCutVal(col, conn);
 
 if ((minVal<0) && (maxVal>0))  /* double color: like microarray */
     rg->maxDev = max(fabs(minVal), maxVal);
 else                           /* single color */
     rg->maxDev = maxVal - minVal;
    
 int i = 0;
 struct slName *sl;
 double val;
 for (sl = gh->sampleList; sl; sl = sl->next)
     {
     struct slName * id = slNameNew(getId(conn, labTable, key, sl->name, value));
     char *cellVal = col->cellVal(col, id, conn);
     if (!cellVal)
 	continue;
 
     val = atof(cellVal);
     if (minCutVal)
 	if (val < atof(minCutVal))
 	    continue;
     if (maxCutVal)
 	if (val > atof(maxCutVal))
 	    continue;
     absVal = fabs(val);
     absVal = reverseColor * (absVal + offset);
 
     if (!(minVal < 0 && maxVal > 0))  
 	absVal = absVal - minVal;
 
     if (val < 0)
 	absVal = -1.0 * absVal;
 
     /* store final data for drawing in ring */
     struct slDouble *sd = slDoubleNew(absVal);
 
     hashAdd(rg->data, sl->name, sd);
     i++;  // increment expId pointer.
 
     freez(&id); 
     }
 }
 
 void ringSetDefaultColor(struct ring *rg)
 {
 if (!rg)
     return;
 
 struct rgbColor *low = NULL;
 struct rgbColor *zero = NULL;
 struct rgbColor *high = NULL;
 
 // default, expression data color scheme:
 if (sameString(rg->type, "CNV") || sameString(rg->type, "SNP"))
     { // CNV, SNP data
     AllocVar(zero);  // white
     zero->r = 255;
     zero->g = 255;
     zero->b = 255;
 
     AllocVar(high);  // red
     high->r = 255;
     high->g = 0;
     high->b = 0;
 
     AllocVar(low);  // blue
     low->r = 0;
     low->g = 0;
     low->b = 255;
     }
 else if (sameString(rg->type, "default"))
     { // default feature color
     AllocVar(zero);  // black
     zero->r = 0;
     zero->g = 0;
     zero->b = 0;
 
     AllocVar(high);  // yellow
     high->r = 255;
     high->g = 255;
     high->b = 0;
 
     AllocVar(low);  // green
     low->r = 0;
     low->g = 255;
     low->b = 0;
     }
 else
     { // expression data, others...
     AllocVar(zero);  // black
     zero->r = 0;
     zero->g = 0;
     zero->b = 0;
 
     AllocVar(high);  // red
     high->r = 255;
     high->g = 0;
     high->b = 0;
 
     AllocVar(low);  // green
     low->r = 0;
     low->g = 255;
     low->b = 0;
     }
 
 rg->zeroColor = zero;
 rg->highColor = high;
 rg->lowColor = low;
 } 
 
 Color getColor(double val, double gain, double maxDev, 
 	       Color *upShades, Color *downShades)
 {
 double absVal = fabs(val);
 double colorScale = COLOR_SCALE / maxDev;
 absVal = absVal * gain;
 
 int colorIndex = (int)(absVal * (EXPR_DATA_SHADES-1.0) * colorScale);
 /* Clip color index to fit inside of array, since we may have brightened it. */
 if (colorIndex < 0)
     colorIndex = 0;
 if (colorIndex >= EXPR_DATA_SHADES)
     colorIndex = EXPR_DATA_SHADES-1;
 
 if(val >= 0)
     return upShades[colorIndex];
 else
     return downShades[colorIndex];
 }
 
 void drawCenter(struct vGfx *vg, struct circleMap *cm)
 {
 int x0 = cm->width / 2;
 int y0 = cm->width / 2;
 Color valCol;
 double tStart, tStop;
 
 Color upShades[EXPR_DATA_SHADES];
 Color downShades[EXPR_DATA_SHADES];    
 
 double numSamples = (double) slCount(cm->samples);
 
 int rStart = 0;
 int rStop = cm->minR;
 
 cm->data = hashNew(0);
 
 struct ring *rg;
 for (rg = cm->rings; rg; rg = rg->next)
     {    
     ringSetDefaultColor(rg);
     
     vgMakeColorGradient(vg, rg->zeroColor, rg->highColor, EXPR_DATA_SHADES, upShades);
     vgMakeColorGradient(vg, rg->zeroColor, rg->lowColor, EXPR_DATA_SHADES, downShades);
 
     double i = 0.0;
     struct slName *sl;
     for (sl = cm->samples; sl; sl = sl->next)
 	{
 	tStart = 2.0 * PI * i / numSamples;
 	tStop  = 2.0 * PI * (i + 1.0) / numSamples;
 	i += 1.0;
 	
 	struct hashEl *el = hashLookup(rg->data, sl->name);
 	if (!el)
 	    continue;  // data exists for sample in ring
 	
 	el = hashLookup(cm->data, sl->name);
 	if (el)
 	    continue;  // only save the first
 	
 	valCol = getColor(rg->median, rg->gain, rg->maxDev, upShades, downShades);
 	vgRadiusBox(vg, x0, y0, rStart, rStop, tStart, tStop, valCol);
 	hashAddInt(cm->data, sl->name, 1);
 	}
     }
 }
 
 void drawRing(struct vGfx *vg, struct circleMap *cm, struct ring *rg)
 {
 if (!rg->data) // nothing to draw
     return;
 
 Color valCol;
 
 Color upShades[EXPR_DATA_SHADES];
 Color downShades[EXPR_DATA_SHADES];
 
 ringSetDefaultColor(rg);
 
 vgMakeColorGradient(vg, rg->zeroColor, rg->highColor, EXPR_DATA_SHADES, upShades);
 vgMakeColorGradient(vg, rg->zeroColor, rg->lowColor, EXPR_DATA_SHADES, downShades);
 
 int x0 = cm->width / 2;
 int y0 = cm->width / 2;
 
 double rStart, rStop, rStep = (rg->maxR - rg->minR) / (double) rg->numElements;
 double tStart, tStop;  // thetas
 
 double numSamples = (double) slCount(cm->samples);
 
 double i = 0.0;
 struct slName *sl;
 struct slDouble *sd, *sdList;
-Color lightGray = vgFindColorIx(vg, 225, 225, 225); 
+
+Color missingColor;
+if (cm->subgrouped)
+    missingColor = MG_WHITE;  // when subgrouped, gray is distracting
+else
+    missingColor = vgFindColorIx(vg, 225, 225, 225);  // light gray
+
+
 for (sl = cm->samples; sl; sl = sl->next)
     {
     tStart = 2.0 * PI * i / numSamples;
     tStop  = 2.0 * PI * (i + 1.0) / numSamples;
     i += 1.0;
     rStart = rg->minR;
     rStop = rg->minR + rStep;
 
     struct hashEl *el = hashLookup(rg->data, sl->name);
     if (!el)
 	{
-	vgRadiusBox(vg, x0, y0, rg->minR, rg->maxR, tStart, tStop, lightGray);
+	vgRadiusBox(vg, x0, y0, rg->minR, rg->maxR, tStart, tStop, missingColor);
 	continue;
 	}
     sdList = el->val;  // slDouble list of expScores
 
     if (rg->drawSummary)
 	{
 	valCol = getColor(rg->median, rg->gain, rg->maxDev, upShades, downShades);
 	vgRadiusBox(vg, x0, y0, rg->minR, rg->maxR, tStart, tStop, valCol);
 	}
     else
 	{
 	for (sd = sdList; sd; sd = sd->next)
 	    {
 	    valCol = getColor(sd->val, rg->gain, rg->maxDev, upShades, downShades);
 	    vgRadiusBox(vg, x0, y0, rStart, rStop, tStart, tStop, valCol);
 	    
 	    rStart += rStep;
 	    rStop += rStep;
 	    }
 	}
     }
 }
 
 void drawRings(struct vGfx *vg, struct circleMap *cm)
 {
 struct ring *rg;
 for (rg = cm->rings; rg; rg = rg->next)
     drawRing(vg, cm, rg);
 }
 
 void drawCircleMap(struct vGfx *vg, struct circleMap *cm)
 {
 if (!cm)
     return;
 
 if (cm->drawSummary)
     drawCenter(vg, cm);
 
 drawRings(vg, cm);
 /* Draw gene label in center */
 MgFont *font = mgSmallFont(); 
 if (mgFontStringWidth(font, cm->name) < cm->minR * 2.0)
     vgTextCentered(vg, 0, 0, cm->width, cm->width, MG_BLACK, font, cm->name);  
 }
 
 char *circleMapGif(struct circleMap *cm)
 /* Create genome GIF file and HT that includes it. */
 {
 if (!cm)
     return NULL;
 
 if (cm->width == 0)
     return NULL;
 
 struct tempName md5Tn;
 char modifier[128];
 safef(modifier, sizeof(modifier), "circleMapGif");
 char *strToHash = cartSettingsString(hgh2Prefix, modifier);
 
 trashDirMD5File(&md5Tn, "hgh", ".gif", strToHash);
 
 off_t size = fileSize(md5Tn.forCgi);
 if (!fileExists(md5Tn.forCgi) || (size == 0) || DEBUG_IMG)
     {
     struct hvGfx *vg = hvGfxOpenGif(cm->width, cm->width, md5Tn.forCgi);
     drawCircleMap(vg->vg, cm);
     hvGfxClose(&vg);
     }
 
 char *filename = replaceChars(md5Tn.forHtml, "..", "");
 return filename;
 }
 
 void layoutCircleMap(struct circleMap *cm, int width)
 {
 if (!cm)
     return;
 cm->width = width;
 cm->maxR = (double) cm->width / 2.0;
 cm->minR = 0.75 * cm->maxR;
 
-int numRings = slCount(cm->rings);
+struct ring *rg;
+int numRings = 0;
+for (rg = cm->rings; rg; rg = rg->next)
+    if (rg->visible)
+	numRings++;
+
 double step = (cm->maxR - cm->minR) / (double) numRings;
 
+double buffer = 0.0;
+if (cm->subgrouped)
+    buffer = step / 5.0;  // put buffer between subgrouped rings
+buffer = 0.0;
+
 double rStart = cm->minR;
-double rStop = cm->minR + step;
+double rStop = cm->minR + (step - buffer);
 
-struct ring *rg;
 for (rg = cm->rings; rg; rg = rg->next)
     {
+    if (!rg->visible)  // keep minR/maxR at zero for no display
+	continue;
     rg->minR = rStart;
     rg->maxR = rStop;
 
     rStart += step;
-    rStop += step;
+    rStop = rStart + (step - buffer);
     }
 }
 
 void setupSummary(struct circleMap *cm)
 {
 char *summarize = cartOptionalString(cart, hgh2Summarize);
 if (!summarize)
     return;
 
 char *summaryDisplay = cartUsualString(cart, hgh2SummaryDisplay, "ring");
 struct slName *names = slNameListFromComma(summarize);
 
 if (sameString(summaryDisplay, "center"))
     cm->drawSummary = TRUE;
 else
     {
     struct ring *rg;
     for (rg = cm->rings; rg; rg = rg->next)
 	{
 	if (!slNameInList(names, rg->name))
 	    continue;
 	rg->drawSummary = TRUE;
 	}
     }
 }
 
 int getSortDirection(char *sortDir)
 {
 if (!sortDir)
     return 1;
 if (sameString(sortDir, "ASC"))
     return 1;
 if (sameString(sortDir, "DESC"))
     return -1;  
 return 0;
 }
 
 void setupGeneSort(struct circleMap *cm)
 {
 char *sortGene = cartOptionalString(cart, hgh2SortGene);
 char *sortDir = cartOptionalString(cart, hgh2SortDir); 
 if (!sortGene)
     return;
 
 struct ring *rg;
 for (rg = cm->rings; rg; rg = rg->next)
     {
     if (!sameString(rg->name, sortGene))
 	continue;
     rg->sort = TRUE;
     rg->sortDir = getSortDirection(sortDir);    
     }
 }
 
 void setupFeatureSort(struct circleMap *cm)
 {
 char *featureSort = cartOptionalString(cart, hgh2FeatureSort);
 char *featureSortDir = cartOptionalString(cart, hgh2FeatureSortDir);
 
 struct slName *s, *sorts = NULL;
 struct slName *d, *sortDirs = NULL;
 if (!featureSort || !featureSortDir)
     return;
 
 sorts = slNameListFromComma(featureSort);
 sortDirs = slNameListFromComma(featureSortDir); 
 if (slCount(sorts) != slCount(sortDirs))
     {
     sorts = NULL;
     sortDirs = NULL;
     }
 
 if (!sorts || !sortDirs)
     return;
 
 struct ring *rg;
 for (rg = cm->rings; rg; rg = rg->next)
     {
     for (s = sorts, d = sortDirs; s && d; s = s->next, d = d->next)
 	{
 	if (!sameString(rg->name, s->name))
 	    continue;
 	rg->sort = TRUE;
 	rg->sortDir = getSortDirection(d->name);
 	}
     }
 }
 
 void sortCircleMap(struct circleMap *cm)
 { /* Sorting goes from inner ring out, averaging all probes within a ring for value */
+if (!cm->samples)
+    return;
 
 /* Setup sorting by gene value */
 setupGeneSort(cm);
 
 /* Setup sorting by features */
 setupFeatureSort(cm);
 
 struct slName *sl;
 struct sortNode *sn, *snList = NULL;
 struct sortList *snl;
 struct ring *rg;
 double val = 0.0;
 for (sl = cm->samples; sl; sl = sl->next)
     {
     AllocVar(sn);
     sn->name = cloneString(sl->name);
     sn->list = NULL;
     for (rg = cm->rings; rg; rg = rg->next)
 	{
-	if (!rg->sort)
+	if (!rg->sort || !rg->data)
 	    continue;
 	struct hashEl *el = hashLookup(rg->data, sn->name);
 	if (!el)
 	    val = rg->maxVal + 0.1;  // if no val, set all to above inner ring's max
 	else
 	    {
 	    struct slDouble *sd, *sdList = el->val;
 	    val = 0.0;
 	    for (sd = sdList; sd; sd = sd->next)
 		val += sd->val;
 	    val = val / (double) slCount(sdList);
 	    }
 	AllocVar(snl);
 	snl->val = val;
 	snl->sortDirection = rg->sortDir;
 	slAddHead(&sn->list, snl);
 	}
     slReverse(&sn->list);
     slAddHead(&snList, sn);
     }
 slSort(&snList, sortNodeCmp);
-
 if (slCount(snList) != slCount(cm->samples))
     errAbort("Sorting failed");
 
 /* Reset sample order for circleMap according to sorted list of nodes */
 slNameFreeList(&cm->samples);
 cm->samples = NULL;
 for (sn = snList; sn; sn = sn->next)
     slNameAddHead(&cm->samples, sn->name);
 slReverse(&cm->samples);
 }
 
 void addGeneToCircleMap(struct circleMap *cm, char *datasets, char *gene)
 {
 if (!ghHash || !gene) // global, make sure it exists
     return;
 
 /* make a geneset containing single gene, for getChromHeatmapHash(...) */
 struct geneSet *gs = AllocA(struct geneSet);
 gs->name = NULL;
 gs->genes = slNameNew(gene);
 
 struct slName *sl, *slList = slNameListFromComma(datasets);
 for (sl = slList; sl; sl = sl->next)
     {
     struct hashEl *el = hashLookup(ghHash, sl->name);
     if (!el)
 	continue;
 
     struct genoHeatmap *gh = el->val;
 
     struct hash *geneHash = NULL;
     getChromHeatmapHash(&geneHash, gh->database, gh->probeTable,
 			gh->name, NULL, gs);    
     if (!geneHash)   // gene hash could not be allocated, bad.
 	continue;
     
-    struct ring *rg = addRingToCircleMap(cm, gene, gh);
+    struct ring *rg = addRingToCircleMap(cm, gene, gh, TRUE);
     addGeneDataToRing(rg, geneHash, gene, gh->sampleList); 
 
     hashFree(&geneHash);
     }
 slReverse(&cm->rings);
 }
 
-void addFeatureToCircleMap(struct circleMap *cm, char *datasets)
+void addFeaturesToCircleMap(struct circleMap *cm, char *datasets)
 {
 char *featureList = cartOptionalString(cart, hgh2FeatureList); 
-if (!featureList)
+char *featureSort = cartOptionalString(cart, hgh2FeatureSort); 
+
+if (!featureList && !featureSort)
     return;
 
+boolean isVisible = TRUE;
+if (!featureList)
+    { /* in order for sorting to work, need features in rings,
+       * but should be invisible unless featureList was set */
+    isVisible = FALSE;
+    featureList = featureSort;
+    }
+
 struct genoHeatmap *gh = NULL;
 struct slName *sl, *slList = slNameListFromComma(datasets);
 for (sl = slList; sl; sl = sl->next)
     {
     struct hashEl *el = hashLookup(ghHash, sl->name);
     if (!el)
 	continue;
     gh = el->val;
     break;
     }
 
 if (!gh)
     return;
 
 struct sqlConnection *conn = hAllocConnProfile(heatMapDbProfile, gh->patDb); 
 if (!conn)
     return;
 
 struct slName *f, *features = slNameListFromComma(featureList);
 /* necessary to reverse so that features show up in order requested
  * can't do a reverse of the rings in case there's already gene rings */
 slReverse(&features);
 
 char *raName = gh->raFile;  
 struct column *col, *colList = getColumns(conn, raName, gh->patDb); 
 for (f = features; f; f = f->next)
     {
     for (col = colList; col != NULL; col = col->next)
 	if (sameString(col->name, f->name))
 	    break;
     if (!col)
 	continue;
     
-    struct ring *rg = addRingToCircleMap(cm, col->name, gh); 
+    struct ring *rg = addRingToCircleMap(cm, col->name, gh, isVisible); 
     addFeatureDataToRing(conn, rg, col, gh);
     }
 hFreeConn(&conn);
 }
 
 void addSamplesToCircleMap(struct circleMap *cm, struct slName *slList)
 {
 if (!cm || !slList)
     return;
 
 struct slName *sl;
 for (sl = slList; sl; sl = sl->next)
     {
     if (slNameInList(cm->samples, sl->name))
 	continue;
     slNameAddHead(&cm->samples, sl->name);
     }
 slReverse(&cm->samples);
 }
 
 void initSamplesCircleMap(struct circleMap *cm, char *datasets)
 {
 if (!ghHash) // global, make sure it exists
     return;
 
 struct slName *sl, *slList = slNameListFromComma(datasets);
 for (sl = slList; sl; sl = sl->next)
     {
     struct hashEl *el = hashLookup(ghHash, sl->name);
     if (!el)
 	continue;
 
     struct genoHeatmap *gh = el->val;
     defaultOrder(gh);   // populate sampleList, ensuring original order (no sorting)
     addSamplesToCircleMap(cm, gh->sampleList);
     }
 }
 
 boolean needSplit(struct ring *rg, struct slName **ptSubsets, int subsetNum)
 {
+if (!rg->data)   // no data, no split!
+    return FALSE;
+
 int subset;
 struct slName *sl;
 
 int needSplit = 0;
 for (subset = 0; subset < subsetNum; subset++)
     {
     for (sl = ptSubsets[subset]; sl; sl = sl->next)
 	{
 	struct hashEl *el = hashLookup(rg->data, sl->name);
 	if (!el)
 	    continue;
 	needSplit++;
+	break;
 	}
     }
 
 return (needSplit > 1);
 }
 
+void removeSamplesFromRing(struct circleMap *cm, struct ring *rg)
+{
+if (!rg->data)  // nothing to remove
+    return;
+
+struct slName *sl, *slList = NULL;
+for (sl = cm->samples; sl; sl = sl->next)
+    {
+    struct hashEl *el = hashLookup(rg->data, sl->name);
+    if (!el)  // not in ring, so keep sample name in circle's list 
+	slNameAddHead(&slList, sl->name);
+    }
+slReverse(&slList);            // maintain original order
+
+slNameFreeList(&cm->samples);  // remove old sample list
+cm->samples = slList;          // set to new list.
+}
+
+
+boolean removeRing(struct circleMap *cm, struct ring *toRemove)
+{
+/* Delete any circleMap samples on the to-be-removed ring */
+removeSamplesFromRing(cm, toRemove);
+
+if (cm->rings == toRemove) /* ring to remove is first, set to next */
+    {
+    cm->rings = toRemove->next;
+    return TRUE;
+    }
+
+struct ring *rg, *nextRg;
+for (rg = cm->rings; rg; rg = rg->next)
+    {
+    nextRg = rg->next;
+    if (nextRg == toRemove)
+	{
+	rg->next = nextRg->next;
+	nextRg->next = NULL;
+	return TRUE;
+	}
+    }
+
+return FALSE;
+}
+
 void splitRing(struct circleMap *cm, struct ring *rg, 
 	       struct slName **ptSubsets, int subsetNum)
 {
-return;
+int subset;
+struct ring *rg2, *ptRings[subsetNum];
+for (subset = 0; subset < subsetNum; subset++)
+    ptRings[subset] = cloneRing(rg);
+
+struct slName *sl;
+for (subset = 0; subset < subsetNum; subset++)
+    {
+    rg2 = ptRings[subset];
+    rg2->data = hashNew(0);
+
+    for (sl = ptSubsets[subset]; sl; sl = sl->next)
+	{
+	struct hashEl *el = hashLookup(rg->data, sl->name);
+	if (!el)
+	    continue;
+	struct slDouble *sdList = el->val;
+
+	hashAdd(rg2->data, sl->name, sdList);
+	hashRemove(rg->data, sl->name);
+	}
+    }
+
+/* Remove empty ring that was split */
+boolean removed = removeRing(cm, rg);
+if (!removed)
+    errAbort("Problem with removing ring");
+
+/* Add new rings */
+for (subset = 0; subset < subsetNum; subset++)
+    slAddHead(&cm->rings, ptRings[subset]);
 }
 
 void splitCircleMap(struct circleMap *cm, struct slName **ptSubsets, int subsetNum)
 {
-struct ring *rg;
-
-for (rg = cm->rings; rg; rg = rg->next)
+struct ring *rg, *nextRg;
+for (rg = cm->rings; rg; rg = nextRg)
     {
+    nextRg = rg->next;
     if (!needSplit(rg, ptSubsets, subsetNum))
 	continue;
     splitRing(cm, rg, ptSubsets, subsetNum);
     }    
 }
 
 
 void setupSubgroups(struct circleMap *cm, char *datasets)
 {
 if (!ghHash) // global, make sure it exists
     return;
 int subsetNum = cartUsualInt(cart, hgh2SubgroupNum, hgh2SubgroupDefaultNum);
 
 struct slName *sl, *slList = slNameListFromComma(datasets);
 for (sl = slList; sl; sl = sl->next)
     {
     struct hashEl *el = hashLookup(ghHash, sl->name);
     if (!el)
 	continue;
     struct genoHeatmap *gh = el->val;
 
     setSubgroups(gh);
 
     char *raName = gh->raFile;
     /* This is the key function to call for subgrouping */
     struct slName **ptSubsets = NULL;
     int ifSubsets = getSubsetsIfAny(gh, subsetNum, raName, &ptSubsets);
     
     if (!ifSubsets)
 	continue;
 
     splitCircleMap(cm, ptSubsets, subsetNum);
+    cm->subgrouped = TRUE;
     }
 }
 
 
 void setupCircleMap(struct circleMap *cm, char *datasets, char *gene, int width)
 {
 /* Add samples for all datasets, list will feature each sample id exactly once*/
 initSamplesCircleMap(cm, datasets);
 
 /* Add gene data (if gene is set) */
 addGeneToCircleMap(cm, datasets, gene);
 
 /* Add feature data (if featureList is set) */
-addFeatureToCircleMap(cm, datasets);
+addFeaturesToCircleMap(cm, datasets);
 
 /* Set up any subgrouping info */
 setupSubgroups(cm, datasets);
 
-/* Layout circle map, setting radial widths */
-layoutCircleMap(cm, width);
-
 /* Setup any summary (displaying only median value) */
 setupSummary(cm);
 
 /* Sort circle map data */
 sortCircleMap(cm);
+
+/* Layout circle map, setting radial widths */
+layoutCircleMap(cm, width);
 }
 
 void modeGetCircleMap()
 {
 char *datasets = cartOptionalString(cart, hgh2Dataset);
 int width = cartUsualInt(cart, hgh2Width, hgHeatmapDefaultCircleWidth);
 char *gene = cartOptionalString(cart, hgh2Gene);
 char *featureList = cartOptionalString(cart, hgh2FeatureList);
 
 if (!ghHash || !datasets)
     return;
 
 if (!featureList && !gene)  // need to specify something to draw
     return;
 
 char *firstFeature = NULL;
 struct slName *f, *features;
 if (featureList)
     {
     features = slNameListFromComma(featureList);
     f = features;
     firstFeature = cloneString(f->name);
     }
 
 char *name;
 if (!gene)
     name = firstFeature;
 else
     name = gene;
 
 struct circleMap *cm = newCircleMap(name);
 setupCircleMap(cm, datasets, gene, width);
 
 struct json *js = newJson();
 char *gGif = circleMapGif(cm);
 jsonAddString(js, "heatmap", gGif);
 
 hPrintf("%s", js->print(js));
 }