Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/Makefile
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/Makefile	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/Makefile	(revision 32055)
@@ -62,4 +62,6 @@
 $(SRC)/high_speed_objects.$(ARCH).o  \
 $(SRC)/high_speed_utils.$(ARCH).o  \
+$(SRC)/relastro_merge_source.$(ARCH).o  \
+$(SRC)/resort_catalog.$(ARCH).o  \
 $(SRC)/relastroVisual.$(ARCH).o
 
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/include/relastro.h
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/include/relastro.h	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/include/relastro.h	(revision 32055)
@@ -15,5 +15,5 @@
 typedef enum {FIT_NONE, FIT_AVERAGE, FIT_PM_ONLY, FIT_PAR_ONLY, FIT_PM_AND_PAR} FitMode;
 
-typedef enum {TARGET_NONE, TARGET_OBJECTS, TARGET_SIMPLE, TARGET_CHIPS, TARGET_MOSAICS, TARGET_HIGH_SPEED} FitTarget;
+typedef enum {TARGET_NONE, TARGET_OBJECTS, TARGET_SIMPLE, TARGET_CHIPS, TARGET_MOSAICS, TARGET_HIGH_SPEED, TARGET_MERGE_SOURCE} FitTarget;
 
 typedef struct {
@@ -94,4 +94,9 @@
 char   SKY_TABLE[256];
 int    SKY_DEPTH;  /** XXX EAM : depth of catalog tables, fix usage */
+
+unsigned int OBJ_ID_SRC;
+unsigned int CAT_ID_SRC;
+unsigned int OBJ_ID_DST;
+unsigned int CAT_ID_DST;
 
 double SIGMA_LIM;
@@ -210,10 +215,10 @@
 void          getfullregion       PROTO((Image *image, off_t Nimage, GSCRegion *fullregion));
 Image        *getimage            PROTO((off_t N));
-Image        *getimages           PROTO((off_t *N));
+Image        *getimages           PROTO((off_t *N, off_t **line_number));
 void          global_stats        PROTO((Catalog *catalog, int Ncatalog));
 void          initGrid            PROTO((int dX, int dY));
 void          initGridBins        PROTO((Catalog *catalog, int Ncatalog));
 void          initImageBins       PROTO((Catalog *catalog, int Ncatalog, int FULLINIT));
-void          initImages          PROTO((Image *input, off_t N));
+void          initImages          PROTO((Image *input, off_t *line_number, off_t N));
 void          initMosaicBins      PROTO((Catalog *catalog, int Ncatalog));
 void          initMosaicGrid      PROTO((Image *image, off_t Nimage));
@@ -372,2 +377,5 @@
 int applyConstraintsB(Catalog *catalog, off_t i);
 void setupAreaSelection(SkyRegion *region);
+
+int relastro_merge_source ();
+void resort_catalog (Catalog *catalog);
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/CoordOps.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/CoordOps.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/CoordOps.c	(revision 32055)
@@ -13,5 +13,5 @@
   Image *images;
 
-  images = getimages (&N);
+  images = getimages (&N, NULL);
 
   NoldCoords = N;
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/FitPM.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/FitPM.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/FitPM.c	(revision 32055)
@@ -73,8 +73,8 @@
   chisq = 0.0;
   for (i = 0; i < Npts; i++) {
-      Xf = fit[0].Ro + fit[0].uR*T[i];
-      Yf = fit[0].Do + fit[0].uD*T[i];
-      chisq += SQ(X[i] - Xf) / SQ(dX[i]);
-      chisq += SQ(Y[i] - Yf) / SQ(dY[i]);
+    Xf = fit[0].Ro + fit[0].uR*T[i];
+    Yf = fit[0].Do + fit[0].uD*T[i];
+    chisq += SQ(X[i] - Xf) / SQ(dX[i]);
+    chisq += SQ(Y[i] - Yf) / SQ(dY[i]);
   }
   fit[0].Nfit = Npts;
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/FitPMandPar.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/FitPMandPar.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/FitPMandPar.c	(revision 32055)
@@ -9,4 +9,5 @@
   double wx, wy, Wx, Wy, Tx, Ty, Tx2, Ty2, Xs, Ys, XT, YT;
   double PR, PD, PRT, PDT, PRX, PDY, PR2, PD2;
+  double chisq, Xf, Yf;
 
   A = array_init (5, 5);
@@ -95,4 +96,16 @@
   /* get the chisq from the matrix values */
 
+  // add up the chi square for the fit
+  chisq = 0.0;
+  for (i = 0; i < Npts; i++) {
+    Xf = fit[0].Ro + fit[0].uR*T[i] + fit[0].dp*pR[i];
+    Yf = fit[0].Do + fit[0].uD*T[i] + fit[0].dp*pD[i];
+    chisq += SQ(X[i] - Xf) / SQ(dX[i]);
+    chisq += SQ(Y[i] - Yf) / SQ(dY[i]);
+  }
+  fit[0].Nfit = Npts;
+
+  // the reduced chisq is divided by (Ndof = 2*Npts - 5)
+  fit[0].chisq = chisq / (2.0*Npts - 5.0);
   return (TRUE);
 }
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/FixProblemImages.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/FixProblemImages.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/FixProblemImages.c	(revision 32055)
@@ -15,5 +15,5 @@
   ALLOCATE (sublist.filename, char *, 1);
 
-  image = getimages (&Nimage);
+  image = getimages (&Nimage, NULL);
 
   Nbad = 0;
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/ImageOps.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/ImageOps.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/ImageOps.c	(revision 32055)
@@ -11,4 +11,5 @@
 static Image        *image;   // list of available images
 static off_t        Nimage;   // number of available images
+static off_t       *LineNumber; // match of subset to full image table
 
 static int         *Ncatlist;  // catalogs associated with each image
@@ -27,6 +28,8 @@
 # endif
 
-Image *getimages (off_t *N) {
+Image *getimages (off_t *N, off_t **line_number) {
+
   *N = Nimage;
+  if (line_number) *line_number = LineNumber;
   return (image);
 }
@@ -42,9 +45,10 @@
 }
 
-void initImages (Image *input, off_t N) {
+void initImages (Image *input, off_t *line_number, off_t N) {
 
   off_t i;
 
   image = input;
+  LineNumber = line_number;
   Nimage = N;
 
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/StarMaps.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/StarMaps.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/StarMaps.c	(revision 32055)
@@ -40,5 +40,5 @@
   off_t i, Nimages;
 
-  images = getimages(&Nimages);
+  images = getimages(&Nimages, NULL);
 
   ALLOCATE (starmap, StarMap, Nimages);
@@ -64,5 +64,5 @@
   gettimeofday (&start, (void *) NULL);
 
-  images = getimages(&Nimages);
+  images = getimages(&Nimages, NULL);
 
   for (i = 0; i < catalog[0].Nmeasure; i++) {
@@ -90,5 +90,5 @@
   int ix, iy;
 
-  images = getimages(&Nimages);
+  images = getimages(&Nimages, NULL);
 
   for (i = 0; i < Nimages; i++) {
@@ -129,5 +129,5 @@
   double L, M, dLmax, dMmax;
 
-  images = getimages(&Nimages);
+  images = getimages(&Nimages, NULL);
 
   dLmax = dMmax = 0.0;
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateChips.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateChips.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateChips.c	(revision 32055)
@@ -18,5 +18,5 @@
   Nskip = Nmosaic = NnewFit = NoldFit = 0;
 
-  image = getimages (&Nimage);
+  image = getimages (&Nimage, NULL);
 
   // save fit results for summary plot
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateMeasures.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateMeasures.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateMeasures.c	(revision 32055)
@@ -13,5 +13,5 @@
     ID_IMAGE_ASTROM_FEW;
 
-  image = getimages (&Nimage);
+  image = getimages (&Nimage, NULL);
 
   for (i = 0; i < Nimage; i++) {
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateObjects.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateObjects.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateObjects.c	(revision 32055)
@@ -196,5 +196,5 @@
       }
       
-      XVERB = FALSE && (catalog[i].measure[m].dM < 0.01) && (N == 6) && (mode == FIT_PM_ONLY);
+      XVERB = (catalog[i].measure[m].dM < 0.01) && (N == 6) && (mode == FIT_PM_ONLY);
 
       // to judge the quality of the PM and PAR fits, we need to fit all three models and compare Chisq
@@ -225,6 +225,5 @@
 
       if (mode == FIT_PM_AND_PAR) {
-	fprintf (stderr, "parallax fitting is still untested (%s, %d)", __FILE__, __LINE__);
-	exit (2);
+	// fprintf (stderr, "parallax fitting is still untested (%s, %d)\n", __FILE__, __LINE__);
 
 	for (k = 0; k < N; k++) {
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateSimple.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateSimple.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/UpdateSimple.c	(revision 32055)
@@ -8,5 +8,5 @@
   StarData *raw, *ref;
 
-  image = getimages (&Nimage);
+  image = getimages (&Nimage, NULL);
 
   for (i = 0; i < Nimage; i++) {
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/args.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/args.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/args.c	(revision 32055)
@@ -1,4 +1,6 @@
 # include "relastro.h"
 void usage (void);
+void usage_merge_source (void);
+void usage_merge_source_id (char *name);
 
 int args (int argc, char **argv) {
@@ -6,8 +8,33 @@
   int N;
   double trange;
+  char *endptr;
 
   /* possible operations */
   FIT_TARGET = TARGET_NONE;
   FIT_MODE = FIT_AVERAGE;
+
+  if ((N = get_argument (argc, argv, "-merge-source"))) {
+    if (N > argc - 6) usage_merge_source();
+    if (strcmp(argv[N+3], "into")) usage_merge_source();
+    FIT_TARGET = TARGET_MERGE_SOURCE;
+    remove_argument (N, &argc, argv);
+    OBJ_ID_SRC = strtol(argv[N], &endptr, 0);
+    if (*endptr) usage_merge_source_id (argv[N]); 
+    remove_argument (N, &argc, argv);
+    CAT_ID_SRC = strtol(argv[N], &endptr, 0);
+    if (*endptr) usage_merge_source_id (argv[N]); 
+    remove_argument (N, &argc, argv);
+    remove_argument (N, &argc, argv); // remove the 'into'
+    OBJ_ID_DST = strtol(argv[N], &endptr, 0);
+    if (*endptr) usage_merge_source_id (argv[N]); 
+    remove_argument (N, &argc, argv);
+    CAT_ID_DST = strtol(argv[N], &endptr, 0);
+    if (*endptr) usage_merge_source_id (argv[N]); 
+    remove_argument (N, &argc, argv);
+
+    if (argc != 1) usage ();
+    return TRUE;
+  }
+
   if ((N = get_argument (argc, argv, "-update-objects"))) {
     remove_argument (N, &argc, argv);
@@ -307,4 +334,6 @@
   fprintf (stderr, "ERROR: USAGE: relastro -region RA RA DEC DEC\n");
   fprintf (stderr, "       OR:    relastro -catalog (ra) (dec)\n\n");
+  fprintf (stderr, "       OR:    relastro -high-speed [options]\n\n");
+  fprintf (stderr, "       OR:    relastro -merge-source [options]\n\n");
   fprintf (stderr, "  specify one of the following modes: \n");
   fprintf (stderr, "  -update-objects\n");
@@ -340,2 +369,12 @@
 } 
 
+void usage_merge_source_id (char *name) {
+
+  fprintf (stderr, "ERROR: invalid ID %s (remember to prefix 0x to hex IDs)\n", name);
+  exit (2);
+}
+
+void usage_merge_source () {
+  fprintf (stderr, "ERROR: USAGE: relastro -merge-source (objID) (catID) into (objID) (catID)\n");
+  exit (2);
+}
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/initialize.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/initialize.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/initialize.c	(revision 32055)
@@ -10,4 +10,6 @@
   ConfigInit (&argc, argv);
   args (argc, argv);
+
+  if (FIT_TARGET == TARGET_MERGE_SOURCE) return;
 
   /* build a list of accepted photcodes. these will be used by bcatalog to accept or
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/load_images.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/load_images.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/load_images.c	(revision 32055)
@@ -44,14 +44,5 @@
   MARKTIME("  select images: %f sec\n", dtime);
 
-  // generate db->vtable from db->ftable based on the selection
-  // XXX does this simply duplicate the memory needlessly?  we recreate these lines
-  // in reload_images.  If we had saved the line numbers, we could avoid this
-  // vtable points *another* copy of the subset rows
-  // (the later call to 'reload_images' copies the subset elements back on top of 
-  // the rows of the vtable)
-  // gfits_vtable_from_ftable (&db[0].ftable, &db[0].vtable, LineNumber, Nsubset);
-  // MARKTIME("converted ftable to vtable: %f sec\n", dtime);
-
-  initImages (subset, Nsubset);
+  initImages (subset, LineNumber, Nsubset);
   MARKTIME("  init images: %f sec\n", dtime);
 
@@ -68,9 +59,10 @@
 
   Image     *image;
-  off_t     Nimage, Nx, i;
+  off_t     Nimage, Nx, i, *LineNumber;
   VTable    *vtable;
 
-  image = getimages (&Nimage);
+  image = getimages (&Nimage, &LineNumber);
 
+  gfits_vtable_from_ftable (&db[0].ftable, &db[0].vtable, LineNumber, Nimage);
   vtable = &db[0].vtable;
 
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/relastro.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/relastro.c	(revision 32054)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/relastro.c	(revision 32055)
@@ -23,4 +23,10 @@
   if (FIT_TARGET == TARGET_HIGH_SPEED) {
     high_speed_catalogs ();
+    exit (0);
+  }
+
+  /* the object analysis is a separate process iterating over catalogs */
+  if (FIT_TARGET == TARGET_MERGE_SOURCE) {
+    relastro_merge_source ();
     exit (0);
   }
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/relastro_merge_source.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/relastro_merge_source.c	(revision 32055)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/relastro_merge_source.c	(revision 32055)
@@ -0,0 +1,138 @@
+# include "relastro.h"
+
+int relastro_merge_source () {
+
+  int i, m;
+  SkyTable *sky = NULL;
+  Catalog catalog_src, catalog_dst;
+
+  // relastro -merge-source (objID) (catID) into (objID) (catID) 
+  // merge the detections of the first source with those of the second source
+  
+  // tell libdvo the CATDIR
+  dvo_set_catdir(CATDIR);
+
+  // load the current sky table (layout of all SkyRegions) 
+  sky = SkyTableLoadOptimal (CATDIR, NULL, NULL, FALSE, SKY_DEPTH, VERBOSE);
+  SkyTableSetFilenames (sky, CATDIR, "cpt");
+  
+  // we are merging the detections of the src into the dst : OBJ_ID_SRC,CAT_ID_SRC -> OBJ_ID_DST,CAT_ID_DST
+
+  // find the catalog containing the src object:
+  dvo_catalog_init (&catalog_src, TRUE);
+  dvo_catalog_init (&catalog_dst, TRUE);
+  SkyRegion *region_src = NULL;
+  SkyRegion *region_dst = NULL;
+
+  // load data from each region file, only use bright stars
+  for (i = 0; i < sky[0].Nregions; i++) {
+    if (sky[0].regions[i].index == CAT_ID_SRC) {
+      catalog_src.filename = sky[0].filename[i];
+      region_src = &sky[0].regions[i];
+    }
+    if (sky[0].regions[i].index == CAT_ID_DST) {
+      catalog_dst.filename = sky[0].filename[i];
+      region_dst = &sky[0].regions[i];
+    }
+  }    
+
+  if (!catalog_src.filename) {
+    fprintf (stderr, "ERROR: cannot find catalog file matching source catID %x\n", CAT_ID_SRC);
+    exit (2);
+  }
+  if (!catalog_dst.filename) {
+    fprintf (stderr, "ERROR: cannot find catalog file matching destination catID %x\n", CAT_ID_DST);
+    exit (2);
+  }
+
+  // load the source object catalog:
+  catalog_src.catformat = dvo_catalog_catformat (CATFORMAT);    // set the default catformat from config data
+  catalog_src.catmode   = dvo_catalog_catmode (CATMODE);        // set the default catmode from config data
+  catalog_src.catflags  = LOAD_AVES | LOAD_MEAS | LOAD_MISS | LOAD_SECF;
+  catalog_src.Nsecfilt  = GetPhotcodeNsecfilt ();
+
+  if (!dvo_catalog_open (&catalog_src, region_src, VERBOSE, "w")) {
+    fprintf (stderr, "ERROR: failure reading catalog %s\n", catalog_src.filename);
+    exit (1);
+  }
+  if (!catalog_src.Naves_disk) {
+    fprintf (stderr, "no data in %s, error in cat ID?\n", catalog_src.filename);
+    exit (1);
+  }
+
+  // find the source object (by objID)
+  int index_src = -1;
+  for (i = 0; i < catalog_src.Naverage; i++) {
+    if (OBJ_ID_SRC != catalog_src.average[i].objID) continue;
+    index_src = i;
+    break;
+  }
+  if (index_src < 0) {
+    fprintf (stderr, "cannot find source object ID %x\n", OBJ_ID_SRC);
+    exit (1);
+  }
+
+  // save the measures from this object 
+  // Measure *measures_src;
+  // ALLOCATE (measures_src, Measure, catalog_src.average[index_src].Nmeasure);
+  // m = catalog_src.average[index_src].measureOffset;
+  // for (i = 0; i < catalog_src.average[index_src].Nmeasure; i++, m++) {
+  //   memcpy (&measures_src[i], &catalog_src.measure[m], sizeof(Measure));
+  // }
+    
+  if (CAT_ID_SRC == CAT_ID_DST) {
+
+    // find the target objects 
+    int index_dst = -1;
+    for (i = 0; i < catalog_src.Naverage; i++) {
+      if (OBJ_ID_DST != catalog_src.average[i].objID) continue;
+      index_dst = i;
+      break;
+    }
+    if (index_dst < 0) {
+      fprintf (stderr, "cannot find source object ID %x\n", OBJ_ID_DST);
+      exit (1);
+    }
+
+    // repoint the src measures at this object
+    Measure *measures_src;
+    ALLOCATE (measures_src, Measure, catalog_src.average[index_src].Nmeasure);
+    m = catalog_src.average[index_src].measureOffset;
+    for (i = 0; i < catalog_src.average[index_src].Nmeasure; i++, m++) {
+
+      // update objID & catID to match the new source
+      catalog_src.measure[m].objID  = catalog_src.average[index_dst].objID;
+      catalog_src.measure[m].catID  = catalog_src.average[index_dst].catID;
+      catalog_src.measure[m].averef = index_dst;
+
+      // get the instantaneous positions:
+      double R = catalog_src.average[index_src].R - catalog_src.measure[m].dR / 3600.0;
+      double D = catalog_src.average[index_src].D - catalog_src.measure[m].dD / 3600.0;
+
+      // update the offset coordinates to match the new source
+      catalog_src.measure[m].dR = 3600.0*(catalog_src.average[index_dst].R - R);
+      catalog_src.measure[m].dD = 3600.0*(catalog_src.average[index_dst].D - D);
+    }
+
+    // update the count? (updated by resort?)
+    catalog_src.average[index_dst].Nmeasure += catalog_src.average[index_src].Nmeasure;
+
+    // for the moment, don't delete this object...
+    catalog_src.average[index_src].Nmeasure = 0;
+    catalog_src.average[index_src].measureOffset = 0;
+
+    catalog_src.sorted = FALSE;
+
+    // resort the measure table and save
+    resort_catalog (&catalog_src);
+
+    dvo_catalog_save (&catalog_src, VERBOSE);
+    dvo_catalog_unlock (&catalog_src);
+    dvo_catalog_free (&catalog_src);
+
+    exit (0);
+  }
+
+  return TRUE;
+}
+
Index: /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/resort_catalog.c
===================================================================
--- /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/resort_catalog.c	(revision 32055)
+++ /branches/eam_branches/ipp-20110710/Ohana/src/relastro/src/resort_catalog.c	(revision 32055)
@@ -0,0 +1,149 @@
+# include "relastro.h"
+
+# define myAbort(MSG) { fprintf (stderr, "%s\n", MSG); abort(); }
+# define myAssert(LOGIC,MSG) { if (!(LOGIC)) { fprintf (stderr, "%s\n", MSG); abort(); } }
+
+// sort the measure Sequence based on the average Sequence entries
+void SortAveMeasMatch (off_t *MEAS, off_t *AVE, off_t N) {
+
+# define SWAPFUNC(A,B){ off_t tmp_meas; off_t tmp_ave;	\
+    tmp_meas = MEAS[A]; MEAS[A] = MEAS[B]; MEAS[B] = tmp_meas;		\
+    tmp_ave  = AVE[A];  AVE[A]  = AVE[B];  AVE[B]  = tmp_ave;		\
+  }
+# define COMPARE(A,B)(AVE[A] < AVE[B])
+  OHANA_SORT (N, COMPARE, SWAPFUNC);
+# undef SWAPFUNC
+# undef COMPARE
+}
+
+// XXX : where is the time going?  perhaps the ALLOCATE?
+// XXX : I don't thnk his is getting the right answer yet.
+
+void resort_catalog (Catalog *catalog) {
+
+  off_t Naverage, Nmeasure;
+  Measure *measure;
+  Average *average;
+  off_t i, j, N, currentAve;
+
+  off_t *measureSeq = NULL;
+  off_t *averageSeq = NULL;
+  Measure *measureTMP = NULL;
+
+  if (catalog[0].sorted == TRUE) return;
+
+  // struct timeval start, stop;
+  // gettimeofday (&start, NULL);
+
+  /* internal counters */
+  Nmeasure = catalog[0].Nmeasure;
+  Naverage = catalog[0].Naverage;
+
+  measure = catalog[0].measure;
+  average = catalog[0].average;
+  
+  // we have a table of average objects and an unsorted table of measurements.  each measurement
+  // has a reference to the average object sequence (as well as an ID)
+  // measure[i].averef -> average[averef]
+  // measure[i].objID = average[averef].objID
+  // measure[i].catID = average[averef].catID
+
+  // we want a sorted measure array with all averef entries in sequence
+
+  ALLOCATE (measureSeq, off_t,   Nmeasure);
+  ALLOCATE (averageSeq, off_t,   Nmeasure);
+
+  for (i = 0; i < Nmeasure; i++) {
+    measureSeq[i] = i;
+    averageSeq[i] = measure[i].averef;
+    
+    if (catalog[0].catformat >= DVO_FORMAT_PS1_V1) {
+      // earlier formats did not carry the objID or catID, so they are not available (we could assign on load, but we don't)
+      myAssert(average[averageSeq[i]].objID == measure[measureSeq[i]].objID, "object / detection mismatch");
+      myAssert(average[averageSeq[i]].catID == measure[measureSeq[i]].catID, "object / detection mismatch");
+    }
+  }
+  
+  // check that averageSeq is now in order
+  // for (i = 1; i < Nmeasure; i++) {
+  //   if (averageSeq[i] < averageSeq[i-1]) {
+  //     fprintf (stderr, "%d ", (int) i);
+  //   }
+  // }
+  // fprintf (stderr, "\n");
+
+  SortAveMeasMatch(measureSeq, averageSeq, Nmeasure);
+  // MARKTIME("sort : %f sec\n", dtime);
+
+  // check that averageSeq is now in order
+  // for (i = 1; i < Nmeasure; i++) {
+  //   if (averageSeq[i] < averageSeq[i-1]) {
+  //     fprintf (stderr, "%d ", (int) i);
+  //   }
+  // }
+  // fprintf (stderr, "\n");
+
+  // copy the measurements in the sorted order
+  ALLOCATE (measureTMP, Measure, Nmeasure);
+  for (i = 0; i < Nmeasure; i++) {
+    j = measureSeq[i];
+    measureTMP[i] = measure[j];
+  }
+  // MARKTIME("assign measure : %f sec\n", dtime);
+
+  // update the values of average.measureOffset and average.Nmeasure
+  FREE(measure);
+  catalog[0].measure = measureTMP;
+
+  N = 0;
+  currentAve = averageSeq[0];
+  average[currentAve].measureOffset = 0;
+  for (i = 0; i < Nmeasure; i++) {
+    if (averageSeq[i] != currentAve) {
+      // we have hit the next entry in the list
+      average[currentAve].Nmeasure = N;
+      N = 0;
+      currentAve = averageSeq[i];
+      average[currentAve].measureOffset = i;
+    }
+    N++;
+  }
+  // N++;
+  average[currentAve].Nmeasure = N;
+  // MARKTIME("update Nmeasure : %f sec\n", dtime);
+
+  int NmeasureTotal = 0;
+  int measureOffsetOK = TRUE;
+  for (i = 0; i < catalog[0].Naverage; i++) {
+    NmeasureTotal += catalog[0].average[i].Nmeasure;
+    if (VERBOSE && !(NmeasureTotal <= catalog[0].Nmeasure)) {
+      fprintf (stderr, "too few measurements: %d %d %d\n", (int) i, NmeasureTotal, (int) catalog[0].Nmeasure);
+    }
+    measureOffsetOK &= (catalog[0].average[i].measureOffset < catalog[0].Nmeasure);
+    if (VERBOSE && !(catalog[0].average[i].measureOffset < catalog[0].Nmeasure)) {
+      fprintf (stderr, "offset too large: %d %d %d\n", (int) i, catalog[0].average[i].Nmeasure, (int) catalog[0].Nmeasure);
+    }
+    measureOffsetOK &= (catalog[0].average[i].measureOffset + catalog[0].average[i].Nmeasure <= catalog[0].Nmeasure);
+    if (VERBOSE && !(catalog[0].average[i].measureOffset + catalog[0].average[i].Nmeasure <= catalog[0].Nmeasure)) {
+      fprintf (stderr, "orrset + Nmeasure too large: %d + %d > %d %d\n", (int) i, catalog[0].average[i].measureOffset, catalog[0].average[i].Nmeasure, (int) catalog[0].Nmeasure);
+    }
+  }
+
+  if (!measureOffsetOK) {
+    fprintf (stderr, "ERROR: catalog %s has an invalid measureOffset\n", catalog[0].filename);
+  }
+
+  if (NmeasureTotal != catalog[0].Nmeasure) {
+    fprintf (stderr, "ERROR: catalog %s has an invalid Nmeasure\n", catalog[0].filename);
+  }
+
+  // MARKTIME("  match time %9.4f sec for %7lld measures, %6lld average\n", dtime, (long long) Nmeasure, (long long) Naverage);
+
+  catalog[0].sorted = TRUE;
+
+  FREE (measureSeq);
+  FREE (averageSeq);
+
+  return;
+}
+