Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/Makefile
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/Makefile	(revision 37544)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/Makefile	(revision 37545)
@@ -77,4 +77,6 @@
 $(SRC)/share_mean_pos.$(ARCH).o \
 $(SRC)/share_images_pos.$(ARCH).o \
+$(SRC)/ICRF.$(ARCH).o \
+$(SRC)/FrameCorrection.$(ARCH).o \
 $(SRC)/ImagePosIO.$(ARCH).o \
 $(SRC)/ImageTable.$(ARCH).o \
@@ -128,4 +130,6 @@
 $(SRC)/plot_scatter.$(ARCH).o	     \
 $(SRC)/plotstuff.$(ARCH).o	     \
+$(SRC)/ICRF.$(ARCH).o \
+$(SRC)/FrameCorrection.$(ARCH).o \
 $(SRC)/relastroVisual.$(ARCH).o \
 $(SRC)/syncfile.$(ARCH).o \
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/doc/pv2.txt
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/doc/pv2.txt	(revision 37544)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/doc/pv2.txt	(revision 37545)
@@ -1,2 +1,20 @@
+
+2014.10.28
+
+Going over the sequence of steps in relastro.  I'm concerned that the
+iteration is not modifying the measure.R,D values, preventing any
+actual iterative improvements.
+
+* UpdateObjects : takes measure.R,D to calculate average.R,D and
+  motion terms
+
+* ImageOps : 
+  * getImageRaw : takes measure.X,Y, applies image.coords to get (L,M),
+    		  (P,Q), etc, calculates 
+  * getImageRef : takes average.R,D, applies image.coords to get
+    		  (P,Q), (L,M), etc
+  * FitChip : fit raw.X,Y to ref.L,M, plus residual map (dX,dY)
+  
+  *** modified raw.L,M values are then dropped and not used to update measure.R.D
 
 2014.07.03
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/include/relastro.h
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/include/relastro.h	(revision 37544)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/include/relastro.h	(revision 37545)
@@ -142,4 +142,13 @@
 } FitStats;
 
+typedef struct {
+  double scale;
+  double **Roff;
+  double **Doff;
+  double  *dR;
+  int *Nra;
+  int Ndec;
+} FrameCorrectionType;
+
 /* global variables set in parameter file */
 # define DVO_MAX_PATH 1024
@@ -188,5 +197,6 @@
 
 int    USE_FIXED_PIXCOORDS;
-int    USE_GAL_MODEL;
+int    USE_GALAXY_MODEL;
+int    USE_ICRF_CORRECT;
 
 int    RESET;
@@ -419,4 +429,5 @@
 StarData *getImageRef (Catalog *catalog, int Ncatalog, off_t im, off_t *Nstars, CoordMode mode);
 StarData *getImageRaw (Catalog *catalog, int Ncatalog, off_t im, off_t *Nstars, CoordMode mode);
+int setImageRaw (Catalog *catalog, int Ncatalog, off_t im, StarData *raw, off_t Nraw, CoordMode mode);
 off_t getImageByID (off_t ID);
 
@@ -555,2 +566,14 @@
 int fit_map (AstromOffsetMap *map, StarData *raw, StarData *ref, int Npts);
 
+// ICRF QSOs : these must be marked in the database (flag on average, flag on measure)
+
+void ICRFinit ();
+int ICRFsave (int cat, int ave, int meas);
+int ICRFdata (int n, int *cat, int *ave, int *meas);
+int ICRFmax ();
+
+FrameCorrectionType *FrameCorrectionInit (double scale);
+void FrameCorrectionFree (FrameCorrectionType *frame);
+int FrameCorrection (Catalog *catalog, int Ncatalog);
+int FrameCorrectionFit (Catalog *catalog, int Ncatalog, SHterms *dR, SHterms *dD);
+int FrameCorrectionFromSH (FrameCorrectionType *frame, SHterms *dR, SHterms *dD);
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ConfigInit.c
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ConfigInit.c	(revision 37544)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ConfigInit.c	(revision 37545)
@@ -36,5 +36,6 @@
 
   if (!ScanConfig (config, "USE_FIXED_PIXCOORDS", "%d", 0, &USE_FIXED_PIXCOORDS))  USE_FIXED_PIXCOORDS = FALSE;
-  if (!ScanConfig (config, "USE_GAL_MODEL",       "%d", 0, &USE_GAL_MODEL))        USE_GAL_MODEL = FALSE;
+  if (!ScanConfig (config, "USE_GALAXY_MODEL",    "%d", 0, &USE_GALAXY_MODEL))     USE_GALAXY_MODEL = FALSE;
+  if (!ScanConfig (config, "USE_ICRF_CORRECT",    "%d", 0, &USE_ICRF_CORRECT))     USE_ICRF_CORRECT = FALSE;
 
   // force CATDIR to be absolute (so parallel mode will work)
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/FrameCorrection.c
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/FrameCorrection.c	(revision 37545)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/FrameCorrection.c	(revision 37545)
@@ -0,0 +1,197 @@
+# include "relastro.h"
+
+/* The sequence of calibration is something like this:
+
+ * apply the current image.coords to measure.X,Y to get measure.R,D
+ * use the measure.R,D values to find the mean position average.R,D
+ * fit the measure.X,Y positions to average.R,D to get image.coords
+
+ At some point in this sequence, we want to compare our current average.R,D (for specific
+ objects) to the ICRF values for those objects to get a correction set in R,D.  
+
+ Apply that correction to get a new set of average.R,D values.  The impact of the
+ adjustment will then be to modify the image.coords to compensate for the offset.
+
+*/
+
+FrameCorrectionType *FrameCorrectionInit (double scale) {
+
+  int i;
+
+  FrameCorrectionType *frame = NULL;
+  ALLOCATE (frame, FrameCorrectionType, 1);
+  
+  frame->scale = scale; // degrees / pixel
+  frame->Ndec = 180.0 / scale;
+  if (frame->Ndec * scale < 180.0) frame->Ndec ++;
+
+  // binD = (DEC + 89.0) / scale
+  ALLOCATE (frame->Nra, int, frame->Ndec);
+  ALLOCATE (frame->Roff, double *, frame->Ndec);
+  ALLOCATE (frame->Doff, double *, frame->Ndec);
+
+  // I am going to avoid corrections at |D| = 90 
+
+  // It may be sleazy, but since I'm using SH, I'm going to just correct the region at |D|
+  // > 80 with a 2D grid of corrections to X,Y projections
+
+  for (i = 0; i < frame->Ndec; i++) {
+    double D = i * scale - 89.0;
+    frame->dR[i] = scale / cos(RAD_DEG*D);
+    frame->Nra[i] = 360.0 / frame->dR[i];
+    if (frame->Nra[i] * frame->dR[i] < 360.0) frame->Nra[i] ++;
+
+    ALLOCATE (frame->Roff[i], double, frame->Nra[i]);
+    ALLOCATE (frame->Doff[i], double, frame->Nra[i]);
+  }
+  return frame;
+}
+
+void FrameCorrectionFree (FrameCorrectionType *frame) {
+
+  if (!frame) return;
+
+  int i;
+
+  for (i = 0; i < frame->Ndec; i++) {
+    if (frame->Roff[i]) free (frame->Roff[i]);
+    if (frame->Doff[i]) free (frame->Doff[i]);
+  }
+
+  if (frame->Nra)  free (frame->Nra);
+  if (frame->Roff) free (frame->Roff);
+  if (frame->Doff) free (frame->Doff);
+
+  free (frame);
+  
+  return;
+}
+
+int FrameCorrection (Catalog *catalog, int Ncatalog) {
+
+  if (!USE_ICRF_CORRECT) return TRUE;
+
+  // we have only a few thousand ICRF reference points.  I need to have a lookup table to
+  // go from ICRF object to the catalog entry.  
+
+  // Lmax in recipe
+  int Lmax = 20;
+
+  // dR and dD will carry the fit coefficients in RA and DEC
+  SHterms *dRc = SHtermsInit (Lmax);
+  SHterms *dDc = SHtermsInit (Lmax);
+
+  FrameCorrectionFit (catalog, Ncatalog, dRc, dDc);
+
+  double pltscale = 0.5; // degrees per pixel
+  FrameCorrectionType *frame = FrameCorrectionInit (pltscale);
+  
+  FrameCorrectionFromSH (frame, dRc, dDc);
+
+  // XXX write out an image to represent the correction
+
+  // Now apply the correction to all of the average.R,D values
+
+  int i, j;
+  for (i = 0; i < Ncatalog; i++) {
+    for (j = 0; j < catalog[i].Naverage; j++) {
+      
+      double R = catalog[i].average[j].R;
+      double D = catalog[i].average[j].D;
+
+      // special handling for the polar regions...
+
+      int iD = (D + 89.0) / frame->scale;
+      int iR = R / frame->dR[iD];
+
+      double dR = frame->Roff[iD][iR];
+      double dD = frame->Doff[iD][iR];
+
+      catalog[i].average[j].R -= dR;
+      catalog[i].average[j].D -= dD;
+    }
+  }
+  FrameCorrectionFree(frame);
+  return TRUE;
+}
+
+
+int FrameCorrectionFit (Catalog *catalog, int Ncatalog, SHterms *dRc, SHterms *dDc) {
+
+  myAssert (dRc->lmax == dDc->lmax, "dR and dD must match\n");
+
+  // allocate an SHterms structure to hold the Ylm values 
+  SHterms *SH = SHtermsInit (dRc->lmax);
+
+  int Nicrf = ICRFmax();
+
+  int i, j;
+  for (i = 0; i < Nicrf; i++) {
+
+    int cat, meas, ave;
+    ICRFdata (i, &cat, &meas, &ave);
+
+    // I need to generate a collection of values dR,dD(R,D), where dR = average.R - ICRF.R, etc.  
+    // I will then fit the set of dR,dD values to a set of spherical harmonics (up to Lmax).
+    // I will then generate correction images dR,dD(R,D)  
+    // I will then correct all values average.R,D by the interpolated values from the images
+    // average.R' = average.R - dR(R,d), average.D' = average.D - dD(R,d)
+
+    Average *average = &catalog[cat].average[ave];
+    Measure *measure = &catalog[cat].measure[meas]; // MeasureTiny?
+
+    // record these in arcsec or degree?
+    // correct for cos(D) or not?
+    double dR = average->R - measure->R;
+    double dD = average->D - measure->D;
+
+    // get the value of Ylm at this coordinate
+    SHtermsForRD (SH, average->R, average->D);
+    
+    // accumulate the dot products
+    for (j = 0; j < SH->Nterms; j++) {
+      dRc->Fr[j] += dR * SH->Fr[j];
+      dRc->Fi[j] += dR * SH->Fi[j];
+      dDc->Fr[j] += dD * SH->Fr[j];
+      dDc->Fi[j] += dD * SH->Fi[j];
+    }
+  }    
+
+  for (j = 0; j < SH->Nterms; j++) {
+    dRc->Fr[j] *= 4.0*M_PI/Nicrf;
+    dRc->Fi[j] *= 4.0*M_PI/Nicrf;
+    dDc->Fr[j] *= 4.0*M_PI/Nicrf;
+    dDc->Fi[j] *= 4.0*M_PI/Nicrf;
+  }
+  SHtermsFree (SH);
+  return TRUE;
+}
+
+int FrameCorrectionFromSH (FrameCorrectionType *frame, SHterms *dR, SHterms *dD) {
+
+  myAssert (dR->lmax == dD->lmax, "dR and dD must match\n");
+
+  // allocate an SHterms structure to hold the Ylm values 
+  SHterms *SH = SHtermsInit (dR->lmax);
+
+  int i, j, k;
+  // Now I need to apply the correction coefficients to SH images to generate an image
+  for (i = 0; i < frame->Ndec; i++) {
+    double D = i * frame->scale - 89.0;
+    for (j = 0; j < frame->Nra[i]; j++) {
+      double R = j * frame->dR[i];
+      
+      SHtermsForRD (SH, R, D);
+      
+      frame->Roff[i][j] = 0.0;
+      frame->Doff[i][j] = 0.0;
+
+      for (k = 0; k < SH->Nterms; k++) {
+	frame->Roff[i][j] += dR->Fr[k]*SH->Fr[k] + dR->Fi[k]*SH->Fi[k];
+	frame->Doff[i][j] += dD->Fr[k]*SH->Fr[k] + dD->Fi[k]*SH->Fi[k];
+      }
+    }      
+  }
+  SHtermsFree (SH);
+  return TRUE;
+}
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ICRF.c
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ICRF.c	(revision 37545)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ICRF.c	(revision 37545)
@@ -0,0 +1,47 @@
+# include "relastro.h"
+
+int Nicrf = 0;
+int NICRF = 0;
+int *ICRFtoCatalog = NULL;
+int *ICRFtoAverage = NULL;
+int *ICRFtoMeasure = NULL;
+
+void ICRFinit () {
+
+  Nicrf = 0;
+  NICRF = 500;
+  ALLOCATE (ICRFtoCatalog, int, NICRF);
+  ALLOCATE (ICRFtoAverage, int, NICRF);
+  ALLOCATE (ICRFtoMeasure, int, NICRF);
+}
+
+int ICRFsave (int cat, int ave, int meas) {
+  
+  ICRFtoCatalog[Nicrf] = cat;
+  ICRFtoAverage[Nicrf] = ave;
+  ICRFtoMeasure[Nicrf] = meas;
+
+  Nicrf ++;
+  if (Nicrf == NICRF) {
+    NICRF += 500;
+
+    REALLOCATE (ICRFtoCatalog, int, NICRF);
+    REALLOCATE (ICRFtoAverage, int, NICRF);
+    REALLOCATE (ICRFtoMeasure, int, NICRF);
+  }
+  return TRUE;
+}
+
+int ICRFmax () {
+  return Nicrf;
+}
+
+int ICRFdata (int n, int *cat, int *ave, int *meas) {
+
+  myAssert (n >= 0 && n < Nicrf, "out of range");
+
+  *cat  = ICRFtoCatalog[n];
+  *ave  = ICRFtoAverage[n];
+  *meas = ICRFtoMeasure[n];
+  return TRUE;
+}
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ImageOps.c
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ImageOps.c	(revision 37544)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ImageOps.c	(revision 37545)
@@ -809,5 +809,5 @@
 
     // XXX subtract off dR,dD GAL here
-    if (USE_GAL_MODEL) {
+    if (USE_GALAXY_MODEL) {
       ref[i].R += measure[0].RoffGAL;
       ref[i].D += measure[0].DoffGAL;
@@ -848,4 +848,56 @@
   *Nstars = N_onImage[im];
   return (ref);
+}
+
+// return StarData values for detections in the specified image, converting coordinates from the
+// chip positions: X,Y -> L,M -> P,Q -> R,D.  This function is used by the image fitting steps, for
+// which the detections have already been filtered when they were loaded (bcatalog)
+int setImageRaw (Catalog *catalog, int Ncatalog, off_t im, StarData *raw, off_t Nraw, CoordMode mode) {
+
+  off_t i, m, c;
+
+  Coords *moscoords;
+
+  moscoords = NULL;
+  if (mode == MODE_MOSAIC) {
+    moscoords = image[im].coords.mosaic;
+    myAssert (moscoords, "coords.mosaic not defined for image %s (%d)", image[im].name, (int) im);
+  }
+
+  myAssert (Nraw == N_onImage[im], "impossible!");
+
+  for (i = 0; i < N_onImage[im]; i++) {
+    m = ImageToMeasure[im][i];
+    c = ImageToCatalog[im][i];
+
+    // XXX should I use the raw coords or just measure.X,Y -> R,D?
+
+    switch (mode) {
+      case MODE_SIMPLE:
+        /* note that for a Simple image, L,M = P,Q */
+        XY_to_LM (&raw[i].L, &raw[i].M, raw[i].X, raw[i].Y, &image[im].coords);
+        raw[i].P = raw[i].L;
+        raw[i].Q = raw[i].M;
+        LM_to_RD (&raw[i].R, &raw[i].D, raw[i].P, raw[i].Q, &image[im].coords);
+        break;
+      case MODE_MOSAIC:
+        XY_to_LM (&raw[i].L, &raw[i].M, raw[i].X, raw[i].Y, &image[im].coords);
+        XY_to_LM (&raw[i].P, &raw[i].Q, raw[i].L, raw[i].M, moscoords);
+        LM_to_RD (&raw[i].R, &raw[i].D, raw[i].P, raw[i].Q, moscoords);
+        break;
+      default:
+	fprintf (stderr, "error: invalid mode in getImageRaw");
+	abort ();
+    }
+
+    MeasureTiny *measure = &catalog[c].measureT[m];
+    measure->R = raw[i].R;
+    measure->D = raw[i].D;
+    if (catalog[c].measure) {
+      catalog[c].measure[m].R = raw[i].R;
+      catalog[c].measure[m].D = raw[i].D;
+    }      
+  }
+  return TRUE;
 }
 
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/UpdateChips.c
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/UpdateChips.c	(revision 37544)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/UpdateChips.c	(revision 37545)
@@ -132,4 +132,7 @@
       continue;
     } 
+
+    // apply the modified R,D back to the measures
+    setImageRaw (catalog, Ncatalog, i, raw, Nraw, MODE_MOSAIC);
 
     saveCenter (image, &Ro[i], &Do[i], i);
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/UpdateObjects.c
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/UpdateObjects.c	(revision 37544)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/UpdateObjects.c	(revision 37545)
@@ -209,5 +209,5 @@
     double Di = getMeanD (&measure[k], average, secfilt);
 
-    if (USE_GAL_MODEL) {
+    if (USE_GALAXY_MODEL) {
       Ri -= measure[0].RoffGAL;
       Di -= measure[0].DoffGAL;
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/bcatalog.c
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/bcatalog.c	(revision 37544)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/bcatalog.c	(revision 37545)
@@ -61,9 +61,10 @@
   }
 
+  int Nicrf = 0;
+
   /* exclude stars not in range or with too few measurements */
   for (i = 0; i < catalog[0].Naverage; i++) {
     if (catalog[0].average[i].Nmeasure <= SRC_MEAS_TOOFEW) {
       continue;
-      // XXX can't (easily) count this for a given chip
     }
     
@@ -77,4 +78,49 @@
     if (RESET) {
       subcatalog[0].average[Naverage].flags &= ~averageBits;
+    }
+
+    // if we want to do static correction to ICRF, keep the ICRF QSOs 
+    if (USE_ICRF_CORRECT && (subcatalog[0].average[Naverage].flags & ID_ICRF_QSO)) {
+      Nm = 0;
+      for (j = 0; j < catalog[0].average[i].Nmeasure; j++) {
+	  
+	offset = catalog[0].average[i].measureOffset + j;
+	  
+	CopyMeasureToTiny (&subcatalog[0].measureT[Nmeasure], &catalog[0].measure[offset]);
+	subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_SKIP_ASTROM;
+	subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_NOCAL;
+	subcatalog[0].measureT[Nmeasure].averef   = Naverage;
+	if (RESET) {
+	  subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_POOR_ASTROM;
+	  subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_AREA;
+	}
+	  
+	// we need to save the location of the ICRF QSOs in the database
+	if (subcatalog[0].measureT[Nmeasure].dbFlags & ID_MEAS_ICRF_QSO) {
+	  ICRFsave (i, Naverage, Nmeasure);
+	  Nicrf ++;
+	}
+
+	Nmeasure ++;
+	Nm ++;
+	if (Nmeasure == NMEASURE) {
+	  NMEASURE += 1000;
+	  REALLOCATE (subcatalog[0].measureT, MeasureTiny, NMEASURE);
+	}
+      }
+      subcatalog[0].average[Naverage].Nmeasure = Nm;
+      Naverage ++;
+      if (Naverage == NAVERAGE) {
+	NAVERAGE += 50;
+	REALLOCATE (subcatalog[0].average, Average, NAVERAGE);
+	REALLOCATE (subcatalog[0].secfilt, SecFilt, NAVERAGE*Nsecfilt);
+      }
+      continue;
+    }
+
+    // if we want to correct to the galay model, we need to use objects with distances determined by Green & Schlafly
+    if (USE_GALAXY_MODEL && !catalog[0].average[i].Nstarpar) {
+      // we may want to add some starpar quality filters...
+      continue;
     }
 
@@ -109,21 +155,21 @@
       // exclude bogus 
       if (ExcludeBogus) {
-	  double Ri = getMeanR_Big (&catalog[0].measure[offset], &catalog[0].average[i], &catalog[0].secfilt[i*Nsecfilt]);
-	  double Di = getMeanD_Big (&catalog[0].measure[offset], &catalog[0].average[i], &catalog[0].secfilt[i*Nsecfilt]);
-	  coords.crval1 = catalog[0].average[i].R;
-	  coords.crval2 = catalog[0].average[i].D;
-	  double Xi, Yi;
-	  RD_to_XY (&Xi, &Yi, Ri, Di, &coords);
-	  double radius = hypot(Xi, Yi);
-	  if (radius > ExcludeBogusRadius) {
-	      NskipBogus ++;
-	      if (VERBOSE2) {
-		FILE *foutput = fbogus ? fbogus : stderr;
-		char *date = ohana_sec_to_date(catalog[0].measure[offset].t);
-		fprintf (foutput, "exclude bogus: %10.6f %10.6f : %10.6f %10.6f : %6.2f %6.2f : %6.2f : %5d %s\n", catalog[0].average[i].R, catalog[0].average[i].D, Ri, Di, Xi, Yi, radius, catalog[0].measure[offset].photcode, date);
-		free (date);
-	      }
-	      continue;
+	double Ri = getMeanR_Big (&catalog[0].measure[offset], &catalog[0].average[i], &catalog[0].secfilt[i*Nsecfilt]);
+	double Di = getMeanD_Big (&catalog[0].measure[offset], &catalog[0].average[i], &catalog[0].secfilt[i*Nsecfilt]);
+	coords.crval1 = catalog[0].average[i].R;
+	coords.crval2 = catalog[0].average[i].D;
+	double Xi, Yi;
+	RD_to_XY (&Xi, &Yi, Ri, Di, &coords);
+	double radius = hypot(Xi, Yi);
+	if (radius > ExcludeBogusRadius) {
+	  NskipBogus ++;
+	  if (VERBOSE2) {
+	    FILE *foutput = fbogus ? fbogus : stderr;
+	    char *date = ohana_sec_to_date(catalog[0].measure[offset].t);
+	    fprintf (foutput, "exclude bogus: %10.6f %10.6f : %10.6f %10.6f : %6.2f %6.2f : %6.2f : %5d %s\n", catalog[0].average[i].R, catalog[0].average[i].D, Ri, Di, Xi, Yi, radius, catalog[0].measure[offset].photcode, date);
+	    free (date);
 	  }
+	  continue;
+	}
       }
 
@@ -150,6 +196,6 @@
       subcatalog[0].measureT[Nmeasure].averef   = Naverage;
       if (RESET) {
-        subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_POOR_ASTROM;
-        subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_AREA;
+	subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_POOR_ASTROM;
+	subcatalog[0].measureT[Nmeasure].dbFlags &= ~ID_MEAS_AREA;
       }
 
@@ -157,6 +203,6 @@
       Nm ++;
       if (Nmeasure == NMEASURE) {
-        NMEASURE += 1000;
-        REALLOCATE (subcatalog[0].measureT, MeasureTiny, NMEASURE);
+	NMEASURE += 1000;
+	REALLOCATE (subcatalog[0].measureT, MeasureTiny, NMEASURE);
       }
     }
@@ -179,5 +225,5 @@
   assert (Nsecfilt == catalog[0].Nsecfilt);
 
-  // limit the total number of stars in the catalog
+// limit the total number of stars in the catalog
   if (MaxDensityUse) {
     LimitDensityCatalog_ByNmeasureGrid (subcatalog, catalog);
@@ -194,4 +240,6 @@
   }
 
+  fprintf (stderr, "added %d ICRF QSO\n", Nicrf);
+  
   return (TRUE);
 }
@@ -210,7 +258,7 @@
 
 # define SWAPFUNC(A,B){ int tmpI; off_t tmpT;	\
-  tmpI = X[A]; X[A] = X[B]; X[B] = tmpI; \
-  tmpT = Y[A]; Y[A] = Y[B]; Y[B] = tmpT; \
-}
+    tmpI = X[A]; X[A] = X[B]; X[B] = tmpI;	\
+    tmpT = Y[A]; Y[A] = Y[B]; Y[B] = tmpT;	\
+  }
 # define COMPARE(A,B)(X[A] > X[B])
 
@@ -224,5 +272,5 @@
 /* generate a grid in a locally-projected coordinate system, try to select average entries 
    from each grid cell in decending Nmeasure order.
- */ 
+*/ 
 int LimitDensityCatalog_ByNmeasureGrid (Catalog *subcatalog, Catalog *oldcatalog) {
 
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/load_catalogs.c
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/load_catalogs.c	(revision 37544)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/load_catalogs.c	(revision 37545)
@@ -88,4 +88,5 @@
 
   bcatalog_show_skips();
+  fprintf (stderr, "included %d ICRF QSOs\n", ICRFmax()); 
 
   Nstar = 0;
Index: /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/relastro_images.c
===================================================================
--- /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/relastro_images.c	(revision 37544)
+++ /branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/relastro_images.c	(revision 37545)
@@ -74,6 +74,6 @@
 	// XXX is is really possible that I do not update the measure.R,D on each loop?? 
 	UpdateObjects (catalog, Ncatalog); // calculate <R>,<D>
+	FrameCorrection (catalog, Ncatalog);
 	UpdateChips (catalog, Ncatalog);   // measure.X,Y -> R,D, fit image.coords
-	FrameCorrection (catalog, Ncatalog);
 	MARKTIME("update chips: %f sec\n", dtime);
       }