Index: /trunk/Ohana/src/libdvo/Makefile
===================================================================
--- /trunk/Ohana/src/libdvo/Makefile	(revision 8294)
+++ /trunk/Ohana/src/libdvo/Makefile	(revision 8295)
@@ -105,4 +105,7 @@
 	cp $(INC)/$* $(DESTINC)/
 
+tabletest: install
+	gcc -L$(LLIB) -I$(LINC) -o tabletest $(SRC)/tabletest.c -ldvo -lFITS -lohana -lm
+
 clean:
 	rm -f */*~
Index: /trunk/Ohana/src/libdvo/src/skyregion_gsc.c
===================================================================
--- /trunk/Ohana/src/libdvo/src/skyregion_gsc.c	(revision 8294)
+++ /trunk/Ohana/src/libdvo/src/skyregion_gsc.c	(revision 8295)
@@ -12,7 +12,9 @@
 			   "s0000", "s0730", "s1500", "s2230", "s3000", "s3730", "s4500", "s5230", "s6000", "s6730", "s7500", "s8230", "none"};
 
+void SkyTableSort (SkyTable *table);
+
 SkyTable *SkyTableFromGSC (char *filename, int depth, int VERBOSE) {
 
-  int i, j, Nx, Ny, No, Nr, NR, Nlines;
+  int i, j, Nx, Ny, No, Nb, Nr, NR, Nlines;
   int nx, ny, Ne, Ns, Nbox;
   double RA0, RA1, DEC0, DEC1, dR, dD;
@@ -23,4 +25,5 @@
   SkyTable *skytable;
   SkyRegion *regions;
+  SkyTable bandregions;
   FILE *f;
   
@@ -95,5 +98,134 @@
   regions[No].childE = Nr;
 
-  /* level 2 : copy the data from the GSC Region files */
+  /* level 2 : identify the blocks in the DEC bands */
+  No = 1;
+  for (i = 0; i < 25; i++) {
+    if (i == 12) continue;
+
+    Nlines = 0;
+    for (j = 0; j < i; j++) Nlines += DecLines[j];
+
+    /* find and save all GSC Regions in this band */
+    Nb = 0;
+    ALLOCATE (bandregions.regions, SkyRegion, DecLines[i]);
+    ALLOCATE (bandregions.filename, char *, DecLines[i]);
+    for (j = 0; j < DecLines[i]; j++) {
+      strncpy (temp, &ftable.buffer[(j + Nlines)*48], 48);
+      temp[49] = 0;
+      hstgsc_hms_to_deg (&RA0, &RA1, &DEC0, &DEC1, &temp[7]);
+      if (DEC1 < DEC0) SWAP (DEC0, DEC1);
+
+      /* check that we are in correct Dec band */
+      if (DEC1 < regions[No].Dmin) {
+	fprintf (stderr, "error: table mis-match (1)\n");
+	fprintf (stderr, "line: %s\n", temp);
+	fprintf (stderr, "region %d: %f - %f\n", No, regions[No].Dmin, regions[No].Dmax);
+	return (NULL);
+      }
+      if (DEC0 > regions[No].Dmax) {
+	fprintf (stderr, "error: table mis-match (2)\n");
+	fprintf (stderr, "line: %s\n", temp);
+	fprintf (stderr, "region %d: %f - %f\n", No, regions[No].Dmin, regions[No].Dmax);
+	return (NULL);
+      }
+
+      /* regions with RA0 < 360 have RA1 == 0.0 */
+      if (RA1 < RA0) RA1 += 360.0;
+
+      bandregions.regions[Nb].Dmin = DEC0;
+      bandregions.regions[Nb].Dmax = DEC1;
+      bandregions.regions[Nb].Rmin = RA0;
+      bandregions.regions[Nb].Rmax = RA1;
+      bandregions.filename[Nb] = NULL;
+      Nb ++;
+    }
+    bandregions.Nregions = Nb;
+
+    /* sort the band regions by Rmin */
+    SkyTableSort (&bandregions);
+    
+    { 
+      int k, Nt, Nrange, Ndiv, Nset;
+      int *transitions;
+      float dDec, dDcur;
+
+      /* search for transitions in the region dDec */
+      ALLOCATE (transitions, int, 50);
+      Nt = 0;
+      dDec = bandregions.regions[0].Dmax - bandregions.regions[0].Dmin;
+      dDcur = dDec;
+      transitions[0] = 0;
+      Nt++;
+      DEC0 = bandregions.regions[0].Dmin;
+      DEC1 = bandregions.regions[0].Dmax;
+      for (j = 0; j < Nb; j++) {
+	DEC0 = MIN (bandregions.regions[j].Dmin, DEC0);
+	DEC1 = MAX (bandregions.regions[j].Dmax, DEC1);
+	dDec = bandregions.regions[j].Dmax - bandregions.regions[j].Dmin;
+	if (dDec != dDcur) {
+	  dDcur = dDec;
+	  transitions[Nt] = j;
+	  Nt++;
+	}
+      }
+      transitions[Nt] = Nb;
+      Nt++;
+
+      /* subdivide each transition range */
+      for (j = 1; j < Nt; j++) {
+	Ne = transitions[j];
+	Ns = transitions[j-1];
+	dDec = bandregions.regions[Ns].Dmax - bandregions.regions[Ns].Dmin;
+	fprintf (stderr, "trans: %d - %d : ", Ns, Ne);
+	fprintf (stderr, "%f,%f - %f,%f\n", 
+		 bandregions.regions[Ns].Rmin, bandregions.regions[Ne-1].Rmax,
+		 bandregions.regions[Ns].Dmin, bandregions.regions[Ns].Dmax);
+
+	Nrange = transitions[j] - transitions[j-1];
+	Ndiv = (int) (0.5 + 7.5 / dDec);
+	switch (Ndiv) {
+	  case 2:
+	    Nset = 4*2;
+	    break;
+	  case 3:
+	    Nset = 6*3;
+	    break;
+	  case 4:
+	    Nset = 8*4;
+	    break;
+	  default:
+	    fprintf (stderr, "programming error\n");
+	    exit (2);
+	}
+	for (k = Ns; k < Ne; k += Nset) {
+	  regions[Nr].Rmin = bandregions.regions[k].Rmin;
+	  if (k + Nset < Ne) {
+	    regions[Nr].Rmax = bandregions.regions[k+Nset-1].Rmax;
+	  } else {
+	    regions[Nr].Rmax = bandregions.regions[Ne-1].Rmax;
+	  }	
+	  regions[Nr].Dmin = DEC0;
+	  regions[Nr].Dmax = DEC1;
+	  fprintf (stderr, "new: %f - %f :%f - %f\n", 
+		   regions[Nr].Rmin, regions[Nr].Rmax, 
+		   regions[Nr].Dmin, regions[Nr].Dmax);
+	}
+
+	regions[Nr].index    =  Nr;
+	regions[Nr].depth    =  2;
+	regions[Nr].parent   =  No;
+	regions[Nr].child    =  TRUE;
+	regions[Nr].table    =  (depth == 2);
+	regions[Nr].childS   =  0;
+	regions[Nr].childE   =  0;
+      }
+      regions[No].childS = Nr;
+      regions[No].childE = Nr;
+      No ++;
+    }
+  }
+
+  /* level 3 : copy the data from the GSC Region files */
+  /* XXX fix this: level here is now 3 */
   No = 1;
   for (i = 0; i < 25; i++) {
@@ -152,4 +284,5 @@
 
   /* subdivide level 2 into NDIV boxes */
+  /* XXX level here is now 4 */
   Ns = regions[regions[0].childS].childS;
   Ne = regions[regions[0].childE - 1].childE;
@@ -201,2 +334,50 @@
   return (skytable);
 }
+
+/* sort skytable by Rmin */
+void SkyTableSort (SkyTable *table) {
+
+  int i, j, k, l, ir, N;
+  SkyRegion *regions;
+  SkyRegion tempregion;
+  char **filename;
+  char *tempfile;
+  
+  N = table[0].Nregions;
+  regions = table[0].regions;
+  filename = table[0].filename;
+
+  if (N < 2) return;
+  l = N >> 1;
+  ir = N - 1;
+  for (;;) {
+    if (l > 0) {
+      l--;
+      tempregion = regions[l];
+      tempfile   = filename[l];
+    } else {
+      tempregion   = regions[ir];
+      tempfile     = filename[ir];
+      regions[ir]   = regions[0];
+      filename[ir] = filename[0];
+      if (--ir == 0) {
+	regions[0]   = tempregion;
+	filename[0] = tempfile;
+	return;
+      }
+    }
+    i = l;
+    j = (l << 1) + 1;
+    while (j <= ir) {
+      if (j < ir && regions[j].Rmin < regions[j+1].Rmin) ++j;
+      if (tempregion.Rmin < regions[j].Rmin) {
+	regions[i]   = regions[j];
+	filename[i] = filename[j];
+	j += (i=j) + 1;
+      } 
+      else j = ir + 1;
+    }
+    regions[i]   = tempregion;
+    filename[i] = tempfile;
+  }
+}
Index: /trunk/Ohana/src/libdvo/src/tabletest.c
===================================================================
--- /trunk/Ohana/src/libdvo/src/tabletest.c	(revision 8295)
+++ /trunk/Ohana/src/libdvo/src/tabletest.c	(revision 8295)
@@ -0,0 +1,10 @@
+# include "dvo.h"
+
+int main (int argc, char **argv) {
+
+  SkyTable *table;
+
+  SkyTableFromGSC (argv[1], 2, TRUE);
+
+  exit (0);
+}