Index: /branches/eam_branches/ipp-20150625/Ohana/src/relastro/include/relastro.h
===================================================================
--- /branches/eam_branches/ipp-20150625/Ohana/src/relastro/include/relastro.h	(revision 38654)
+++ /branches/eam_branches/ipp-20150625/Ohana/src/relastro/include/relastro.h	(revision 38655)
@@ -294,4 +294,6 @@
 int    CHIPORDER;
 int    CHIPMAP;
+
+int    N_BOOTSTRAP_SAMPLES;
 
 int MaxDensityUse;
@@ -710,5 +712,5 @@
 int CatalogMaxNmeasure (Catalog *catalog, int Ncatalog);
 int FitAstromPoints_Project (FitStats *fitStats, double *Tmean, double *Trange, double *parRange);
-int UpdateObjects_SelectMeasures (FitStats *fit, Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int cat, off_t measOff, int isStack);
+int UpdateObjects_SelectMeasures (FitStats *fit, Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int isStack);
 int UpdateObjects_Stack (Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int Nsecfilt, FitStats *fitStats);
 int UpdateObjects_Chips (Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int Nsecfilt, FitStats *fitStats, int cat, off_t measOff);
@@ -722,2 +724,3 @@
 void FitStatsReset (FitStats *tgt);
 FitStats *FitStatsInit (int Nmax, int Nboot);
+int FitAstromResultSetPM (FitAstromResult *fit, int Nfit, Average *average);
Index: /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/FitPM.c
===================================================================
--- /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/FitPM.c	(revision 38654)
+++ /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/FitPM.c	(revision 38655)
@@ -1,4 +1,5 @@
 # include "relastro.h"
 
+// initial values of *fit are ignored
 int FitPM (FitAstromResult *fit, FitAstromData *data, FitAstromPoint *points, int Npoints) {
 
Index: /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/FitPMandPar.c
===================================================================
--- /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/FitPMandPar.c	(revision 38654)
+++ /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/FitPMandPar.c	(revision 38655)
@@ -1,4 +1,5 @@
 # include "relastro.h"
 
+// initial values of *fit are ignored
 int FitPMandPar (FitAstromResult *fit, FitAstromData *data, FitAstromPoint *points, int Npoints) {
 
Index: /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/UpdateObjectOffsets.c
===================================================================
--- /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/UpdateObjectOffsets.c	(revision 38654)
+++ /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/UpdateObjectOffsets.c	(revision 38655)
@@ -185,4 +185,6 @@
     // XXX note that the above pass in the flag as decimal -- also note that args.c cannot handle 0xHEX values
 
+    if (N_BOOTSTRAP_SAMPLES > 1) strextend (&command, "-bootstrap-samples %d", N_BOOTSTRAP_SAMPLES); 
+
     if (DCR_BLUE_COLOR_POS && DCR_BLUE_COLOR_NEG) {
       strextend (&command, "-dcr-blue-color %s %s", DCR_BLUE_COLOR_POS, DCR_BLUE_COLOR_NEG); 
Index: /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/UpdateObjects.c
===================================================================
--- /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/UpdateObjects.c	(revision 38654)
+++ /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/UpdateObjects.c	(revision 38655)
@@ -1,5 +1,4 @@
 # include "relastro.h"
 # define PAR_TOOFEW 5
-# define NBOOT 100
 
 // This function operates on both Measure and MeasureTiny.  In the big stages, this should
@@ -19,6 +18,6 @@
   FitStats *sumStatsStack = FitStatsInit (0, 0);
 
-  FitStats *fitStatsChips = FitStatsInit (NmeasureMax, NBOOT);
-  FitStats *fitStatsStack = FitStatsInit (NmeasureMax, NBOOT);
+  FitStats *fitStatsChips = FitStatsInit (NmeasureMax, N_BOOTSTRAP_SAMPLES);
+  FitStats *fitStatsStack = FitStatsInit (NmeasureMax, N_BOOTSTRAP_SAMPLES);
 
   int i;
@@ -39,4 +38,9 @@
       SecFilt *secfilt = &catalog[i].secfilt[j*Nsecfilt];
 
+      int XVERB = FALSE;
+      XVERB |= (average[0].objID == OBJ_ID_SRC) && (average[0].catID == CAT_ID_SRC);
+      XVERB |= (average[0].objID == OBJ_ID_DST) && (average[0].catID == CAT_ID_DST);
+      if ((OBJ_ID_SRC || OBJ_ID_DST) && !XVERB) continue;
+
       UpdateObjects_Stack(average, secfilt, measure, measureBig, Nsecfilt, fitStatsStack);
       UpdateObjects_Chips(average, secfilt, measure, measureBig, Nsecfilt, fitStatsChips, i, m);
@@ -90,5 +94,5 @@
 
   // select the measurements to be used in this analysis
-  UpdateObjects_SelectMeasures (fitStats, average, secfilt, measure, measureBig, cat, measOff, FALSE);
+  UpdateObjects_SelectMeasures (fitStats, average, secfilt, measure, measureBig, FALSE);
 
   // if there are no exposure detections, use the stack position
@@ -123,24 +127,24 @@
     }
 
-    // XXX handle the case with only 1 or few points
-
-    fitStats->Nfit = 0;
-    for (k = 0; k < fitStats->NfitAlloc; k++) {
-      BootstrapResample (fitStats->sample, fitStats->points, fitStats->Npoints);
-      if (!FitPM (&fitStats->fit[k], fitStats->fitdataPM, fitStats->sample, fitStats->Npoints)) continue;
-      fitStats->Nfit ++;
-    }
-    BootstrapRobustStats (&fitPM, fitStats->fit, fitStats->Nfit, FIT_RESULT_RA);
-    BootstrapRobustStats (&fitPM, fitStats->fit, fitStats->Nfit, FIT_RESULT_DEC);
-    BootstrapRobustStats (&fitPM, fitStats->fit, fitStats->Nfit, FIT_RESULT_uR);
-    BootstrapRobustStats (&fitPM, fitStats->fit, fitStats->Nfit, FIT_RESULT_uD);
-
+    if (fitStats->NfitAlloc == 1) {
+      // if N_BOOTSTRAP_SAMPLES = 1, no bootstrap resampling:
+      FitPM (&fitPM, fitStats->fitdataPM, fitStats->points, fitStats->Npoints);
+    } else {
+      fitStats->Nfit = 0;
+      for (k = 0; k < fitStats->NfitAlloc; k++) {
+	BootstrapResample (fitStats->sample, fitStats->points, fitStats->Npoints);
+	if (!FitPM (&fitStats->fit[k], fitStats->fitdataPM, fitStats->sample, fitStats->Npoints)) continue;
+	fitStats->Nfit ++;
+      }
+      BootstrapRobustStats (&fitPM, fitStats->fit, fitStats->Nfit, FIT_RESULT_RA);
+      BootstrapRobustStats (&fitPM, fitStats->fit, fitStats->Nfit, FIT_RESULT_DEC);
+      BootstrapRobustStats (&fitPM, fitStats->fit, fitStats->Nfit, FIT_RESULT_uR);
+      BootstrapRobustStats (&fitPM, fitStats->fit, fitStats->Nfit, FIT_RESULT_uD);
+    }
     FitAstromSetChisq (&fitPM, fitStats->points, fitStats->Npoints, FIT_PM_ONLY);
 
     // project Ro, Do back to RA,DEC
     XY_to_RD (&fitPM.Ro, &fitPM.Do, fitPM.Ro, fitPM.Do, &fitStats->coords);
-    if (XVERB) fprintf (stderr, "project: %f %f : %f %f : %f\n", fitPM.Ro, fitPM.Do, fitPM.uR, fitPM.uD, fitPM.p);
     if (fabs(fitPM.Ro) < 0.01) fprintf (stderr, "watch out for 0,360 boundary\n");
-    // XXX : does this make sense at 0,360 boundary?
 
     fitPM.p  = fitPM.dp  = 0.0;
@@ -168,23 +172,28 @@
     }
 
-    fitStats->Nfit = 0;
-    for (k = 0; k < fitStats->NfitAlloc; k++) {
-      BootstrapResample (fitStats->sample, fitStats->points, fitStats->Npoints);
-      FitPMandPar (&fitStats->fit[k], fitStats->fitdataPar, fitStats->sample, fitStats->Npoints);
-      fitStats->Nfit ++;
-    }
-    BootstrapRobustStats (&fitPar, fitStats->fit, fitStats->Nfit, FIT_RESULT_RA);
-    BootstrapRobustStats (&fitPar, fitStats->fit, fitStats->Nfit, FIT_RESULT_DEC);
-    BootstrapRobustStats (&fitPar, fitStats->fit, fitStats->Nfit, FIT_RESULT_uR);
-    BootstrapRobustStats (&fitPar, fitStats->fit, fitStats->Nfit, FIT_RESULT_uD);
-    BootstrapRobustStats (&fitPar, fitStats->fit, fitStats->Nfit, FIT_RESULT_PLX);
-
+    if (fitStats->NfitAlloc == 1) {
+      // if N_BOOTSTRAP_SAMPLES = 1, no bootstrap resampling:
+      FitPMandPar (&fitPar, fitStats->fitdataPar, fitStats->points, fitStats->Npoints);
+    } else {
+      fitStats->Nfit = 0;
+      for (k = 0; k < fitStats->NfitAlloc; k++) {
+	BootstrapResample (fitStats->sample, fitStats->points, fitStats->Npoints);
+	FitPMandPar (&fitStats->fit[k], fitStats->fitdataPar, fitStats->sample, fitStats->Npoints);
+	fitStats->Nfit ++;
+      }
+      BootstrapRobustStats (&fitPar, fitStats->fit, fitStats->Nfit, FIT_RESULT_RA);
+      BootstrapRobustStats (&fitPar, fitStats->fit, fitStats->Nfit, FIT_RESULT_DEC);
+      BootstrapRobustStats (&fitPar, fitStats->fit, fitStats->Nfit, FIT_RESULT_uR);
+      BootstrapRobustStats (&fitPar, fitStats->fit, fitStats->Nfit, FIT_RESULT_uD);
+      BootstrapRobustStats (&fitPar, fitStats->fit, fitStats->Nfit, FIT_RESULT_PLX);
+    }
     FitAstromSetChisq (&fitPar, fitStats->points, fitStats->Npoints, FIT_PM_AND_PAR);
 
+    // project Ro, Do back to RA,DEC
     XY_to_RD (&fitPar.Ro, &fitPar.Do, fitPar.Ro, fitPar.Do, &fitStats->coords);
+    if (fabs(fitPar.Ro) < 0.01) fprintf (stderr, "watch out for 0,360 boundary\n");
+
     average[0].flags |= ID_STAR_FIT_PAR;
     fitStats->Npar ++;
-
-    if (fabs(fitPar.Ro) < 0.01) fprintf (stderr, "watch out for 0,360 boundary\n");
 
     // XXX a hard-wired hack...
@@ -199,12 +208,18 @@
     // if we only have one point, this is silly...
     
-    fitStats->Nfit = 0;
-    for (k = 0; k < fitStats->NfitAlloc; k++) {
-      BootstrapResample (fitStats->sample, fitStats->points, fitStats->Npoints);
-      FitPosPMfixed (&fitStats->fit[k], fitStats->fitdataPos, fitStats->sample, fitStats->Npoints);
-      fitStats->Nfit ++;
-    }
-    BootstrapRobustStats (&fitPos, fitStats->fit, fitStats->Nfit, FIT_RESULT_RA);
-    BootstrapRobustStats (&fitPos, fitStats->fit, fitStats->Nfit, FIT_RESULT_DEC);
+    if (fitStats->NfitAlloc == 1) {
+      FitAstromResultSetPM (&fitPos, 1, average);
+      FitPosPMfixed (&fitPos, fitStats->fitdataPos, fitStats->points, fitStats->Npoints);
+    } else {
+      fitStats->Nfit = 0;
+      FitAstromResultSetPM (fitStats->fit, fitStats->NfitAlloc, average);
+      for (k = 0; k < fitStats->NfitAlloc; k++) {
+	BootstrapResample (fitStats->sample, fitStats->points, fitStats->Npoints);
+	FitPosPMfixed (&fitStats->fit[k], fitStats->fitdataPos, fitStats->sample, fitStats->Npoints);
+	fitStats->Nfit ++;
+      }
+      BootstrapRobustStats (&fitPos, fitStats->fit, fitStats->Nfit, FIT_RESULT_RA);
+      BootstrapRobustStats (&fitPos, fitStats->fit, fitStats->Nfit, FIT_RESULT_DEC);
+    }
     FitAstromSetChisq (&fitPos, fitStats->points, fitStats->Npoints, FIT_AVERAGE);
 
@@ -215,5 +230,5 @@
   }
 
-  // XXX update the bit flags of which points were used
+  // update the bit flags of which points were used
   for (k = 0; k < fitStats->Npoints; k++) {
     int Nm = fitStats->points[k].measure;
@@ -223,17 +238,41 @@
   }
 
-# if (0)
+  // we can set the star reference-image color only if we have loaded the image data
   int setRefColor = areImagesMatched();
   if (setRefColor) {
+    float *C_blue = NULL;
+    float *C_red = NULL;
+    ALLOCATE (C_blue, float, fitStats->Npoints);
+    ALLOCATE (C_red, float, fitStats->Npoints);
+
+    int NcBlue = 0;
+    int NcRed = 0;
+
+    for (k = 0; k < fitStats->Npoints; k++) {
+      int Nm = fitStats->points[k].measure;
+      float colorBlue = getColorBlue (measOff + Nm, cat);
+      if (!isnan(colorBlue)) {
+	C_blue[NcBlue] = colorBlue;
+	NcBlue++;
+      }
+      float colorRed = getColorRed (measOff + Nm, cat);
+      if (!isnan(colorRed)) {
+	C_red[NcRed] = colorRed;
+	NcRed++;
+      }
+    }
+
     // need to reassign here if isfinite()
     float colorMedian;
-    dsort (C_blue, NcBlue);
+    fsort (C_blue, NcBlue);
     colorMedian = (NcBlue > 0) ? C_blue[(int)(0.5*NcBlue)] : NAN;
     average[0].refColorBlue = colorMedian;
-    dsort (C_red, NcRed);
+    fsort (C_red, NcRed);
     colorMedian = (NcRed > 0) ? C_red[(int)(0.5*NcRed)] : NAN;
     average[0].refColorRed = colorMedian;
-  }
-# endif
+
+    free (C_blue);
+    free (C_red);
+  }
 
   /* choose the result based on the chisq values */
@@ -362,5 +401,5 @@
 
   // select the measurements to be used in this analysis
-  UpdateObjects_SelectMeasures (fitStats, average, secfilt, measure, measureBig, 0, 0, TRUE);
+  UpdateObjects_SelectMeasures (fitStats, average, secfilt, measure, measureBig, TRUE);
 
   // too few measurements for average position (require 2 values)
@@ -421,5 +460,5 @@
 }
 
-int UpdateObjects_SelectMeasures (FitStats *fit, Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int cat, off_t measOff, int isStack) {
+int UpdateObjects_SelectMeasures (FitStats *fit, Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int isStack) {
 
   // I've already allocated fit->points (and fit->sample) with space for fit->NpointsAlloc entries
@@ -512,11 +551,4 @@
 
     points[Npoints].dT = measure[k].dt;
-
-    // XXX 
-    int setRefColor = 0;
-    if (setRefColor) {
-      points[Npoints].C_blue = getColorBlue (measOff+k, cat);
-      points[Npoints].C_red  = getColorRed (measOff+k, cat);
-    }
 
     points[Npoints].measure = k;
@@ -720,4 +752,24 @@
       myAbort ("invalid mode");
   }
+  fit->Nfit = Npoints;
   return (TRUE);
 }
+
+int FitAstromResultSetPM (FitAstromResult *fit, int Nfit, Average *average) {
+
+  int i;
+
+  if (USE_GALAXY_MODEL) {
+    for (i = 0; i < Nfit; i++) {
+      fit->uR = average->uRgal;
+      fit->uD = average->uDgal;
+    }
+  } else {
+    for (i = 0; i < Nfit; i++) {
+      fit->uR = 0.0;
+      fit->uD = 0.0;
+    }
+  }
+
+  return TRUE;
+}
Index: /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/args.c
===================================================================
--- /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/args.c	(revision 38654)
+++ /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/args.c	(revision 38655)
@@ -271,4 +271,15 @@
     remove_argument (N, &argc, argv);
     APPLY_OFFSETS = TRUE;
+  }
+
+  N_BOOTSTRAP_SAMPLES = 1;
+  if ((N = get_argument (argc, argv, "-bootstrap-samples"))) {
+    remove_argument (N, &argc, argv);
+    N_BOOTSTRAP_SAMPLES = atoi (argv[N]);
+    remove_argument (N, &argc, argv);
+    if ((N_BOOTSTRAP_SAMPLES < 1) || ((N_BOOTSTRAP_SAMPLES > 1) && (N_BOOTSTRAP_SAMPLES < 20))) {
+      fprintf (stderr, "-bootstrap-samples must be either 1 (no sampling) or >= 20\n");
+      exit (2);
+    }
   }
 
@@ -697,4 +708,15 @@
     remove_argument (N, &argc, argv);
     MaxDensityUse = TRUE;
+  }
+
+  N_BOOTSTRAP_SAMPLES = 1;
+  if ((N = get_argument (argc, argv, "-bootstrap-samples"))) {
+    remove_argument (N, &argc, argv);
+    N_BOOTSTRAP_SAMPLES = atoi (argv[N]);
+    remove_argument (N, &argc, argv);
+    if ((N_BOOTSTRAP_SAMPLES < 1) || ((N_BOOTSTRAP_SAMPLES > 1) && (N_BOOTSTRAP_SAMPLES < 20))) {
+      fprintf (stderr, "-bootstrap-samples must be either 1 (no sampling) or >= 20\n");
+      exit (2);
+    }
   }
 
Index: /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/high_speed_catalogs.c
===================================================================
--- /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/high_speed_catalogs.c	(revision 38654)
+++ /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/high_speed_catalogs.c	(revision 38655)
@@ -176,6 +176,7 @@
     // XXX note that the above pass in the flag as decimal -- also note that args.c cannot handle 0xHEX values
 
-    if (MinBadQF > 0.0)        { strextend (&command, "-min-bad-psfqf %f", MinBadQF); }
-    if (MaxMeanOffset != 10.0) { strextend (&command, "-max-mean-offset  %f", MaxMeanOffset); }
+    if (MinBadQF > 0.0)          strextend (&command, "-min-bad-psfqf %f", MinBadQF); 
+    if (MaxMeanOffset != 10.0)   strextend (&command, "-max-mean-offset  %f", MaxMeanOffset); 
+    if (N_BOOTSTRAP_SAMPLES > 1) strextend (&command, "-bootstrap-samples %d", N_BOOTSTRAP_SAMPLES); 
 
     if (WHERE_A[0]) { strextend (&command, "-D WHERE_A \"%s\"", WHERE_A); }
Index: /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/hpm_catalogs.c
===================================================================
--- /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/hpm_catalogs.c	(revision 38654)
+++ /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/hpm_catalogs.c	(revision 38655)
@@ -174,6 +174,7 @@
     // XXX note that the above pass in the flag as decimal -- also note that args.c cannot handle 0xHEX values
 
-    if (MinBadQF > 0.0)        { strextend (&command, "-min-bad-psfqf %f", MinBadQF); }
-    if (MaxMeanOffset != 10.0) { strextend (&command, "-max-mean-offset  %f", MaxMeanOffset); }
+    if (MinBadQF > 0.0)          strextend (&command, "-min-bad-psfqf %f", MinBadQF); 
+    if (MaxMeanOffset != 10.0)   strextend (&command, "-max-mean-offset  %f", MaxMeanOffset); 
+    if (N_BOOTSTRAP_SAMPLES > 1) strextend (&command, "-bootstrap-samples %d", N_BOOTSTRAP_SAMPLES); 
 
     if (TimeSelect) { 
Index: /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/relastro_objects.c
===================================================================
--- /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/relastro_objects.c	(revision 38654)
+++ /branches/eam_branches/ipp-20150625/Ohana/src/relastro/src/relastro_objects.c	(revision 38655)
@@ -172,6 +172,7 @@
     // XXX note that the above pass in the flag as decimal -- also note that args.c cannot handle 0xHEX values
 
-    if (MinBadQF > 0.0)        { strextend (&command, "-min-bad-psfqf %f", MinBadQF); }
-    if (MaxMeanOffset != 10.0) { strextend (&command, "-max-mean-offset  %f", MaxMeanOffset); }
+    if (MinBadQF > 0.0)          strextend (&command, "-min-bad-psfqf %f", MinBadQF); 
+    if (MaxMeanOffset != 10.0)   strextend (&command, "-max-mean-offset  %f", MaxMeanOffset); 
+    if (N_BOOTSTRAP_SAMPLES > 1) strextend (&command, "-bootstrap-samples %d", N_BOOTSTRAP_SAMPLES); 
 
     if (TimeSelect) { 
