0817a98144008545e4ae66bb12a23eed676875a8
braney
  Wed Apr 9 17:52:26 2014 -0700
make snakes more compact
diff --git src/hg/hgTracks/snakeTrack.c src/hg/hgTracks/snakeTrack.c
index 5215073..c42e61a 100644
--- src/hg/hgTracks/snakeTrack.c
+++ src/hg/hgTracks/snakeTrack.c
@@ -7,94 +7,102 @@
 #include "linefile.h"
 #include "jksql.h"
 #include "hdb.h"
 #include "hgTracks.h"
 #include "chainBlock.h"
 #include "chainLink.h"
 #include "chainDb.h"
 #include "chainCart.h"
 #include "errCatch.h"
 #include "twoBit.h"
 #include "bigWarn.h"
 #include <pthread.h>
 #include "trackHub.h"
 #include "limits.h"
 #include "snakeUi.h"
+#include "bits.h"
 
 #include "halBlockViz.h"
 
 // this is the number of pixels used by the target self-align bar
 #define DUP_LINE_HEIGHT	4
 
 struct snakeFeature
     {
     struct snakeFeature *next;
     int start, end;			/* Start/end in browser coordinates. */
     int qStart, qEnd;			/* query start/end */
     int level;				/* level in snake */
     int orientation;			/* strand... -1 is '-', 1 is '+' */
     boolean drawn;			/* did we draw this feature? */
     char *qSequence;			/* may have sequence, or NULL */
     char *tSequence;			/* may have sequence, or NULL */
     char *qName;			/* chrom name on other species */
+    unsigned pixX1, pixX2;              /* pixel coordinates within window */
     };
 
 static int snakeFeatureCmpTStart(const void *va, const void *vb)
 /* sort by start position on the target sequence */
 {
 const struct snakeFeature *a = *((struct snakeFeature **)va);
 const struct snakeFeature *b = *((struct snakeFeature **)vb);
 int diff = a->start - b->start;
 
 return diff;
 }
 
 static int snakeFeatureCmpQStart(const void *va, const void *vb)
 /* sort by start position on the query sequence */
 {
 const struct snakeFeature *a = *((struct snakeFeature **)va);
 const struct snakeFeature *b = *((struct snakeFeature **)vb);
 int diff = a->qStart - b->qStart;
 
 if (diff == 0)
     {
     diff = a->start - b->start;
     }
 
 return diff;
 }
 
-
 // static machinery to calculate full and pack snake
 
+#define NUM_LEVELS	1000
+
 struct level
 {
 boolean init;		/* has this level been initialized */
 int orientation;	/* strand.. see above */
 unsigned long edge;	/* the leading edge of this level */
 int adjustLevel;	/* used to compress out the unused levels */
 boolean hasBlock;	/* are there any blocks in this level */
+
+// these fields are used to compact full snakes
+unsigned *connectsToLevel;  /* what levels does this level connect to */
+Bits *pixels;               /* keep track of what pixels are used */
+struct snakeFeature *blocks; /* the ungapped blocks attached to this level */
 };
 
-static struct level Levels[10000]; /* for packing the snake, not re-entrant! */
+static struct level Levels[NUM_LEVELS]; /* for packing the snake, not re-entrant! */
 static int maxLevel = 0;     	     /* deepest level */	
 
 /* blocks that make it through the min size filter */
 static struct snakeFeature *newList = NULL;   
 
 static void clearLevels()
-/* clear out the data structure that we use to calculate full snakes */
+/* clear out the data structure that we use to calculate snakes */
 {
 int ii;
 
 for(ii=0; ii < sizeof(Levels) / sizeof(Levels[0]); ii++)
     Levels[ii].init = FALSE;
 maxLevel = 0;
 }
 
 static void calcFullSnakeHelper(struct snakeFeature *list, int level)
 // calculate a full snake
 // updates global newList with unsorted blocks that pass the min
 // size filter.
 {
 struct snakeFeature *cb = list;
 struct snakeFeature *proposedList = NULL;
@@ -200,33 +208,204 @@
 
 // transfer proposedList to new block list
 for(temp=proposedList; temp; temp = next)
     {
     next = temp->next;
     temp->next = NULL;
     slAddHead(&newList, temp);
     }
 }
 
 struct snakeInfo
 {
 int maxLevel;
 } snakeInfo;
 
+static void freeFullLevels()
+// free the connection levels and bitmaps for each level
+{
+int ii;
+
+for(ii=0; ii <= maxLevel; ii++)
+    {
+    freez(&Levels[ii].connectsToLevel);
+    bitFree(&Levels[ii].pixels);
+    }
+}
+
+static void initializeFullLevels(struct snakeFeature *sfList)
+/* initialize levels for compact snakes */
+{
+int ii;
+
+for(ii=0; ii <= maxLevel; ii++)
+    {
+    Levels[ii].connectsToLevel = needMem((maxLevel+1) * sizeof(unsigned));
+    Levels[ii].pixels = bitAlloc(insideWidth);
+    Levels[ii].blocks = NULL;
+    }
+
+struct snakeFeature *sf, *prev = NULL, *next;
+int prevLevel = -1;
+double scale = scaleForWindow(insideWidth, winStart, winEnd);
+
+for(sf=sfList; sf; prev = sf, sf = next)
+    {
+    next = sf->next;
+
+    int sClp = (sf->start < winStart) ? winStart : sf->start;
+    sf->pixX1 = round((sClp - winStart)*scale);
+    int eClp = (sf->end > winEnd) ? winEnd : sf->end;
+    sf->pixX2 = round((eClp - winStart)*scale);
+
+    bitSetRange(Levels[sf->level].pixels, sf->pixX1, sf->pixX2 - sf->pixX1);
+
+    if (prevLevel != -1)
+	Levels[sf->level].connectsToLevel[prevLevel] = 1;
+
+    if(next != NULL)
+	Levels[sf->level].connectsToLevel[next->level] = 1;
+
+    prevLevel = sf->level;
+
+    slAddHead(&Levels[sf->level].blocks, sf);
+    }
+}
+
+static int highestLevelBelowUs(int level)
+// is there a level below us that doesn't connect to us
+{
+int newLevel;
+
+// find the highest level below us that we connect to
+for(newLevel = level - 1; newLevel >= 0; newLevel--)	 
+    if (Levels[level].connectsToLevel[newLevel])
+	break;
+
+// if the next level is more than one below us, we may
+// be able to push our blocks to it.
+if ((newLevel < 0) || (newLevel + 1 == level))
+    return -1;
+
+return newLevel + 1;
+}
+
+// the number of pixels we need clear on either side of the blocks we push
+#define EXTRASPACE	10
+
+static boolean checkBlocksOnLevel(int level, struct snakeFeature *sfList)
+// is there enough pixels on this level to hold our blocks?
+{
+struct snakeFeature *sf;
+
+for(sf=sfList; sf;  sf = sf->next)
+    {
+    if (bitCountRange(Levels[level].pixels, sf->pixX1 - EXTRASPACE, (sf->pixX2 - sf->pixX1) + 2*EXTRASPACE))
+	return FALSE;
+    }
+return TRUE;
+}
+
+static void collapseLevel(int oldLevel, int newLevel)
+// push level up to higher level
+{
+struct snakeFeature *sf;
+struct snakeFeature *next;
+
+for(sf=Levels[oldLevel].blocks; sf;  sf = next)
+    {
+    next = sf->next;
+    sf->level = newLevel;
+    slAddHead(&Levels[newLevel].blocks, sf);
+
+    bitSetRange(Levels[sf->level].pixels, sf->pixX1, sf->pixX2 - sf->pixX1);
+    }
+Levels[oldLevel].blocks = NULL;
+Levels[oldLevel].pixels = bitAlloc(insideWidth);
+}
+
+static void remapConnections(int oldLevel, int newLevel)
+// if we've moved a level, we need to remap all the connections
+{
+int ii, jj;
+
+for(ii=0; ii <= maxLevel; ii++)
+    {
+    for(jj=0; jj <= maxLevel; jj++)
+	{
+	if (Levels[ii].connectsToLevel[oldLevel])
+	    {
+	    Levels[ii].connectsToLevel[oldLevel] = 0;
+	    Levels[ii].connectsToLevel[newLevel] = 1;
+	    }
+	}
+    }
+}
+
+static struct snakeFeature *compactSnakes(struct snakeFeature *sfList)
+// collapse levels that will fit on a higer level
+{
+int ii;
+
+initializeFullLevels(sfList);
+
+// start with third level to see what can be compacted
+for(ii=2; ii <= maxLevel; ii++) 
+    {
+    int newLevel;
+    // is there a level below us that doesn't connect to us
+    if ((newLevel = highestLevelBelowUs(ii)) > 0)
+	{
+	int jj;
+	for (jj=newLevel; jj < ii; jj++)
+	    {
+	    if (checkBlocksOnLevel(jj, Levels[ii].blocks))
+		{
+		collapseLevel(ii, jj);
+		remapConnections(ii, jj);
+		break;
+		}
+	    }
+	}
+    }
+
+// reattach blocks in the levels to linkedFeature
+struct snakeFeature *newList = NULL;
+for(ii=0; ii <= maxLevel; ii++)
+    {
+    newList = slCat(Levels[ii].blocks, newList);
+    }
+slSort(&newList, snakeFeatureCmpQStart);
+
+// figure out the new max
+int newMax = 0;
+for(ii=0; ii <= maxLevel; ii++)
+    if (Levels[ii].blocks != NULL)
+	newMax = ii;
+
+freeFullLevels();
+maxLevel = newMax;
+
+return newList;
+}
+
 static void calcFullSnake(struct track *tg, void *item)
+// calculate a full snake
 {
 struct linkedFeatures  *lf = (struct linkedFeatures *)item;
+
+// if there aren't any blocks, don't bother
 if (lf->components == NULL)
     return;
 
 // we use the codons field to keep track of whether we already
 // calculated the height of this snake
 if (lf->codons == NULL)
     {
     clearLevels();
     struct snakeFeature *sf;
 
     // this will destroy lf->components, and add to newList
     calcFullSnakeHelper((struct snakeFeature *)lf->components, 0);
     lf->components = (struct simpleFeature *)newList;
     newList = NULL;
     slSort(&lf->components, snakeFeatureCmpQStart);
@@ -274,30 +453,32 @@
     // now figure out how to remap the blocks
     int ii;
     int count = 0;
     for(ii=0; ii < oldMax + 1; ii++)
 	{
 	Levels[ii].adjustLevel = count;
 	if ((Levels[ii].init) && (Levels[ii].hasBlock))
 	    count++;
 	}
     maxLevel = count - 1;
 
     // remap blocks
     for(sf=(struct snakeFeature *)lf->components; sf; sf = sf->next)
 	sf->level = Levels[sf->level].adjustLevel;
 
+    // now compact the snakes
+    lf->components = (void *)compactSnakes((struct snakeFeature *)lf->components);
     struct snakeInfo *si;
     AllocVar(si);
     si->maxLevel = maxLevel;
     lf->codons = (struct simpleFeature *)si;
     }
 }
 
 static void calcPackSnakeHelper(struct snakeFeature *list, int level)
 // calculate a packed snake
 // updates global newList with unsorted blocks that pass the min
 // size filter.
 {
 struct snakeFeature *cb = list;
 struct snakeFeature *proposedList = NULL;