Index: /branches/eam_branches/ohana.20190329/src/addstar/include/ukirt_uhs.h
===================================================================
--- /branches/eam_branches/ohana.20190329/src/addstar/include/ukirt_uhs.h	(revision 40832)
+++ /branches/eam_branches/ohana.20190329/src/addstar/include/ukirt_uhs.h	(revision 40832)
@@ -0,0 +1,50 @@
+
+// measure[0] = g, measure[1] = b, measure[2] = r, 
+
+typedef struct {
+  Average average;
+  Measure measure[3];
+  int flag; // in a subset?
+  int found; // assigned to an object?
+} UKIRT_UHS_Stars;
+
+int   HOST_ID;
+char *HOSTDIR;
+char *CPT_FILE;
+char *INPUT;
+
+AddstarClientOptions args_loadukirt_uhs (int *argc, char **argv, AddstarClientOptions options);
+// AddstarClientOptions args_loadukirt_uhs_client (int *argc, char **argv, AddstarClientOptions options);
+
+int loadukirt_uhs_table (int Nstart, int Nend, SkyList *skylistInput, HostTable *hosts, char **filename, AddstarClientOptions *options);
+
+UKIRT_UHS_Stars *loadukirt_uhs_make_subset (UKIRT_UHS_Stars *stars, int Nstars, int start, SkyRegion *region, int *nsubset);
+
+// int loadukirt_uhs_save_remote (UKIRT_UHS_Stars *stars, int Nstars, HostTable *hosts, SkyRegion *region, char *fullname, AddstarClientOptions *options);
+// int save_remote_host (HostInfo *host);
+
+// int init_remote_hosts (void);
+// void free_remote_hosts (void);
+// int find_empty_slot (void);
+// int harvest_all (void);
+// int harvest_host (void);
+
+int loadukirt_uhs_catalog (UKIRT_UHS_Stars *stars, int Nstars, SkyRegion *region, char *filename, AddstarClientOptions *options);
+
+// int galactic_to_celestial (double *R, double *D, double l, double b);
+
+int find_matches_ukirt_uhs (SkyRegion *region, UKIRT_UHS_Stars *stars, int Nstars, Catalog *catalog, AddstarClientOptions *options);
+
+int loadukirt_uhs_save_stars (char *filename, UKIRT_UHS_Stars *stars, int Nstars);
+// UKIRT_UHS_Stars *loadukirt_uhs_load_stars (char *filename, int *nstars);
+
+UKIRT_UHS_Stars *loadukirt_uhs_readstars (char *filename, UKIRT_UHS_Stars *stars, int *nstars, AddstarClientOptions *options);
+
+int loadukirt_uhs_sortStars (UKIRT_UHS_Stars *stars, int Nstars);
+
+// int loadukirt_uhs_tmpdir (void);
+
+//  roll these into a function with a structure to carry the Nlast, etc?
+char *dparse_csv_rpt (double *X, int Nwant, int Nlast, char *line, int *status);
+char *iparse_csv_rpt (int *X,    int Nwant, int Nlast, char *line, int *status);
+char *jparse_csv_rpt (uint64_t *X, int Nwant, int Nlast, char *line, int *status);
Index: /branches/eam_branches/ohana.20190329/src/addstar/src/args_loadukirt_uhs.c
===================================================================
--- /branches/eam_branches/ohana.20190329/src/addstar/src/args_loadukirt_uhs.c	(revision 40832)
+++ /branches/eam_branches/ohana.20190329/src/addstar/src/args_loadukirt_uhs.c	(revision 40832)
@@ -0,0 +1,137 @@
+# include "addstar.h"
+# include "ukirt_uhs.h"
+
+static void help (void);
+
+AddstarClientOptions args_loadukirt_uhs (int *argc, char **argv, AddstarClientOptions options) {
+  
+  int N;
+
+  /* check for help request */
+  if (get_argument (*argc, argv, "-help") ||
+      get_argument (*argc, argv, "-h")) {
+    help ();
+  }
+
+  /* basic mode: image, list, refcat */
+  options.mode = ADDSTAR_MODE_REFCAT;
+
+  /* we do not allow a subset to be extracted -- all or nothing, babe */
+  UserPatch.Rmin = 0;
+  UserPatch.Rmax= 360;
+  UserPatch.Dmin = -90;
+  UserPatch.Dmax = +90;
+
+  // XXX for the moment, make this selection manual.  it needs to be automatic 
+  // based on the state of the SkyTable
+  HOST_ID = 0;
+  PARALLEL = FALSE;
+  if ((N = get_argument (*argc, argv, "-parallel"))) {
+    PARALLEL = TRUE;
+    remove_argument (N, argc, argv);
+  }
+  // this is a test mode : rather than launching the remote jobs and waiting for completion,
+  // relphot will simply list the remote command and wait for the user to signal completion
+  PARALLEL_MANUAL = FALSE;
+  if ((N = get_argument (*argc, argv, "-parallel-manual"))) {
+    PARALLEL = TRUE; // -parallel-manual implies -parallel
+    PARALLEL_MANUAL = TRUE;
+    remove_argument (N, argc, argv);
+  }
+  // this is a test mode : rather than launching the relphot_client jobs remotely, they are 
+  // run in serial via 'system'
+  PARALLEL_SERIAL = FALSE;
+  if ((N = get_argument (*argc, argv, "-parallel-serial"))) {
+    if (PARALLEL_MANUAL) {
+      fprintf (stderr, "ERROR: cannot mix -parallel-manual and -parallel-serial\n");
+      exit (1);
+    }
+    PARALLEL = TRUE; // -parallel-serial implies -parallel
+    PARALLEL_SERIAL = TRUE;
+    remove_argument (N, argc, argv);
+  }
+
+  /* only add to existing regions */
+  options.existing_regions = FALSE;
+  if ((N = get_argument (*argc, argv, "-existing-regions"))) {
+    options.existing_regions = TRUE;
+    remove_argument (N, argc, argv);
+  }
+  /* only add to existing objects */
+  options.only_match = FALSE;
+  if ((N = get_argument (*argc, argv, "-only-match"))) {
+    options.only_match = TRUE;
+    remove_argument (N, argc, argv);
+  }
+  /* replace measurement, don't duplicate (ref/cat only) */
+  options.replace = FALSE;
+  if ((N = get_argument (*argc, argv, "-replace"))) {
+    options.replace = TRUE;
+    remove_argument (N, argc, argv);
+  }
+  /* override any header PHOTCODE values */
+  options.photcode = 0;
+
+  /* accept proper-motion (& parallax) data from reference */
+  ACCEPT_MOTION = FALSE;
+  if ((N = get_argument (*argc, argv, "-accept-motion"))) {
+    ACCEPT_MOTION = TRUE;
+    remove_argument (N, argc, argv);
+  }
+
+  /* provide a time for dataset */
+  options.timeref = 0; 
+  if ((N = get_argument (*argc, argv, "-time"))) {
+    time_t tmp;
+    remove_argument (N, argc, argv);
+    if (!ohana_str_to_time (argv[N], &tmp)) { 
+      fprintf (stderr, "syntax error in time\n");
+      exit (1);
+    }
+    options.timeref = tmp;
+    remove_argument (N, argc, argv);
+  }
+
+  /* extra error messages */
+  VERBOSE = FALSE;
+  if ((N = get_argument (*argc, argv, "-v"))) {
+    VERBOSE = TRUE;
+    remove_argument (N, argc, argv);
+  }
+
+  /* other addstar options which cannot be used in loadukirt_uhs */
+  // options.timeref = 0; 
+  options.mosaic = FALSE;
+  options.skip_missed = FALSE;
+  options.closest = FALSE;
+  options.nosort = FALSE;
+  options.update = FALSE;
+  options.only_images = FALSE;
+  options.calibrate = FALSE;
+  options.quality_airmass = FALSE;
+  ACCEPT_ASTROM = FALSE;
+  FORCE_READ = FALSE;
+  TEXTMODE = FALSE;
+  SUBPIX = FALSE;
+  DUMP = NULL;
+
+  if (*argc < 2) {
+    fprintf (stderr, "USAGE: loadukirt_uhs [options] (fitsfile) [..more files]\n");
+    exit (2);
+  }
+  return (options);
+}
+
+static void help () {
+
+  fprintf (stderr, "USAGE: loadukirt_uhs [options] (file) [..more files]\n");
+  fprintf (stderr, "  add data from stellar parameter file to DVO\n\n");
+
+  fprintf (stderr, "  optional flags:\n");
+  fprintf (stderr, "  -only-match           	  : only add measurements to existing objects\n");
+  fprintf (stderr, "  -replace              	  : replace time/photcode measurements (no duplication)\n");
+  fprintf (stderr, "  -v                    	  : verbose mode\n");
+  fprintf (stderr, "  -help                 	  : this list\n");
+  fprintf (stderr, "  -h                    	  : this list\n\n");
+  exit (2);
+}
Index: /branches/eam_branches/ohana.20190329/src/addstar/src/find_matches_ukirt_uhs.c
===================================================================
--- /branches/eam_branches/ohana.20190329/src/addstar/src/find_matches_ukirt_uhs.c	(revision 40832)
+++ /branches/eam_branches/ohana.20190329/src/addstar/src/find_matches_ukirt_uhs.c	(revision 40832)
@@ -0,0 +1,346 @@
+# include "addstar.h"
+# include "ukirt_uhs.h"
+
+int find_matches_ukirt_uhs (SkyRegion *region, UKIRT_UHS_Stars *stars, int NstarsIn, Catalog *tgtcat, AddstarClientOptions *options) {
+
+  int Nsecfilt = GetPhotcodeNsecfilt ();
+  int Nsec     = GetPhotcodeNsec (options->photcode);
+
+  /** allocate local arrays (stars) **/
+  ALLOCATE_PTR (X1, double, NstarsIn);
+  ALLOCATE_PTR (Y1, double, NstarsIn);
+  ALLOCATE_PTR (N1, off_t,  NstarsIn);
+
+  /** allocate local arrays (tgtcat) **/
+  off_t NAVE = tgtcat[0].Naverage;
+  off_t Nave = tgtcat[0].Naverage;
+  ALLOCATE_PTR (X2, double, NAVE);
+  ALLOCATE_PTR (Y2, double, NAVE);
+  ALLOCATE_PTR (N2, off_t,  NAVE);
+
+  /* for secfilt j and star i, secfilt[i*Nsecfilt+j] */
+
+  /* internal counters */
+  off_t Nmatch = 0;
+  off_t NMEAS = tgtcat[0].Nmeasure;
+  off_t Nmeas = tgtcat[0].Nmeasure;
+
+  off_t *next_meas = NULL;
+
+  // current max obj ID for this tgtcat
+  unsigned int objID = tgtcat[0].objID;
+  unsigned int catID = tgtcat[0].catID;
+
+  /* project onto rectilinear grid with 1 arcsec pixels. the choice of ARC projection has
+   * the advantage that every point in R,D has a mapping to a unique X,Y.  However, note
+   * that not all possible X,Y points map back to R,D and the local plate scale changes
+   * substantially far from the projection pole. We use the center of the region (tgtcat)
+   * for crval1,2.
+   */
+
+  Coords tcoords;
+  InitCoords (&tcoords, "DEC--ARC");
+  tcoords.crval1 = 0.5*(region[0].Rmin + region[0].Rmax);
+  if ((region[0].Dmax < 90) && (region[0].Dmin > -90)) {
+    tcoords.crval2 = 0.5*(region[0].Dmin + region[0].Dmax);
+  } else {
+    tcoords.crval2 = (region[0].Dmax >= 90) ? 90.0 : -90.0;
+  }
+  tcoords.cdelt1 = tcoords.cdelt2 = 1.0 / 3600.0;
+
+  /* build spatial index (RA sort) referencing input array sequence */
+  off_t Nstars = 0;
+  for (off_t i = 0; i < NstarsIn; i++) {
+    int status = RD_to_XY (&X1[Nstars], &Y1[Nstars], stars[i].average.R, stars[i].average.D, &tcoords);
+    if (!status) continue;
+    N1[Nstars] = i;
+    Nstars ++;
+  }
+  if (Nstars < 1) {
+    if (VERBOSE) fprintf (stderr, "skipping %s, no overlapping stars\n", tgtcat[0].filename);
+    free (X1);
+    free (Y1);
+    free (N1);
+    free (X2);
+    free (Y2);
+    free (N2);
+    return 0;
+  }
+  if (Nstars > 1) sort_coords_index (X1, Y1, N1, Nstars);
+
+  /* build spatial index (RA sort) */
+  for (off_t i = 0; i < Nave; i++) {
+    RD_to_XY (&X2[i], &Y2[i], tgtcat[0].average[i].R, tgtcat[0].average[i].D, &tcoords);
+    N2[i] = i;
+  }
+  if (Nave > 1) sort_coords_index (X2, Y2, N2, Nave);
+
+  /* set up pointers for linked list of measure */
+  if (tgtcat[0].sorted && (tgtcat[0].Nmeasure == tgtcat[0].Nmeasure_disk)) {
+    // this version is only valid if we have done a full tgtcat load, and if the tgtcat
+    // is sorted while processed
+    next_meas = init_measure_links (tgtcat[0].average, Nave, tgtcat[0].measure, Nmeas);
+  } else {
+    next_meas = build_measure_links (tgtcat[0].average, Nave, tgtcat[0].measure, Nmeas);
+  }    
+
+  /* choose a radius for matches */
+  double RADIUS = (options->radius == 0) ? 2.0 : options->radius; /* provided by config */
+  double RADIUS2 = RADIUS*RADIUS;
+
+  /****************** find matched stars ********************/
+
+  for (off_t i = 0, j = 0; (i < Nstars) && (j < Nave); ) {
+    if (!finite(X1[i]) || !finite(Y1[i])) { 
+      i++; 
+      continue;
+    }
+    if (!finite(X2[j]) || !finite(Y2[j])) { 
+      j++; 
+      continue;
+    }
+    
+    /* negative dX: j is too large */
+    double dX = X1[i] - X2[j];
+    if (dX <= -1.02*RADIUS) {
+      i++;
+      continue;
+    }
+    /* positive dX, i is too large */
+    if (dX >= 1.02*RADIUS) {
+      j++;
+      continue;
+    }
+
+    // skip this star if already assigned to an object
+    if (stars[N1[i]].found) {
+        i++;
+        continue;
+    }
+
+    /* this block will match a given detection to the closest object within range of that detection.
+       XXX note that this matches ALL detections within range of the single object to that same object 
+       this is bad, but I cannot just go in linear order (ie, mark off each object as they are
+       used).  I should make a list of all Nobj * Ndet pairs in range and choose the matches
+       based on their separations.  UGH
+     */
+    
+    /* within match range; look for matches */
+    off_t Jmin = -1;
+    double Rmin = RADIUS2;
+    for (off_t J = j; (dX > -1.02*RADIUS) && (J < Nave); J++) {
+      /* find closest match for this detection */
+      dX = X1[i] - X2[J];
+      double dY = Y1[i] - Y2[J];
+      double dR = dX*dX + dY*dY;
+      if (dR > RADIUS2) continue;
+      if (dR < Rmin) {
+	Rmin = dR;
+	Jmin  = J;
+      }
+    }
+
+    /* no match, try next detection */ 
+    if (Jmin == -1) {
+      i++;
+      continue;
+    }
+
+    /*** a match is found, add to average, measure ***/
+    Nmatch ++;
+    off_t n = N2[Jmin];
+    off_t N = N1[i];
+
+    /* make sure there is space for next entry */
+    if (Nmeas >= NMEAS - stars[N].average.Nmeasure) {
+      NMEAS = Nmeas + 1000 + stars[N].average.Nmeasure;
+      REALLOCATE (next_meas, off_t, NMEAS);
+      REALLOCATE (tgtcat[0].measure, Measure, NMEAS);
+    }
+
+    /** add measurements for this star **/
+    for (int mSrc = 0; mSrc < stars[N].average.Nmeasure; mSrc++) {
+
+      /* add to end of measurement list */
+      add_meas_link (&tgtcat[0].average[n], next_meas, Nmeas, NMEAS);
+
+      /** add measurements for this star **/
+
+      // set the new measurements
+      tgtcat[0].measure[Nmeas]          = stars[N].measure[mSrc];
+
+      // measure now carries R,D (not dR,dD) 
+      tgtcat[0].measure[Nmeas].dbFlags  = 0;
+      tgtcat[0].measure[Nmeas].averef   = n;
+      tgtcat[0].measure[Nmeas].objID    = tgtcat[0].average[n].objID;
+      tgtcat[0].measure[Nmeas].catID    = tgtcat[0].catID;
+      
+      float dRoff = dvoOffsetR(&tgtcat[0].measure[Nmeas], &tgtcat[0].average[n]);
+
+      // rationalize dR
+      if (dRoff > +180.0*3600.0) {
+	// average on high end of boundary, move star up
+	tgtcat[0].measure[Nmeas].R += 360.0;
+	dRoff -= 360.0*3600.0;
+      }
+      if (dRoff < -180.0*3600.0) {
+	// average on low end of boundary, move star down
+	tgtcat[0].measure[Nmeas].R -= 360.0;
+	dRoff += 360.0*3600.0;
+      }
+      if (fabs(dRoff) > 10*RADIUS) {
+	// take declination into account and check again.
+	double cosD = cos(RAD_DEG*tgtcat[0].average[n].D);
+	if (fabs(dRoff*cosD) > 10*RADIUS) {
+	  fprintf (stderr, "error: %10.6f,%10.6f vs %10.6f,%10.6f (%f,%f vs %f,%f)\n", 
+		   tgtcat[0].average[n].R, tgtcat[0].average[n].D, 
+		   stars[N].average.R, stars[N].average.D,
+		   X1[i], X2[Jmin], 
+		   Y1[i], Y2[Jmin]);
+	  // XXX abort on this? -- this is a bad failure...
+	}
+      }
+      
+      /* Nm is updated, but not written out in -update mode (for existing entries)
+	 Nm is recalculated in build_meas_links if loaded table is not sorted */
+      tgtcat[0].average[n].Nmeasure ++;
+      Nmeas ++;
+      
+      /* if we choose to flag close encounters, see find_matches.c */
+      /* if we choose to calculate RA,DEC averages, see update_coords.c */
+    }
+    
+    if (ACCEPT_MOTION) {
+      tgtcat[0].average[n].uR    = stars[N].average.uR;
+      tgtcat[0].average[n].uD    = stars[N].average.uD;
+      tgtcat[0].average[n].duR   = stars[N].average.duR;
+      tgtcat[0].average[n].duD   = stars[N].average.duD;
+      tgtcat[0].average[n].P     = stars[N].average.P;
+      tgtcat[0].average[n].dP    = stars[N].average.dP;
+      tgtcat[0].average[n].Tmean = stars[N].average.Tmean;
+    }
+	
+    stars[N].found = TRUE;
+    i++;
+  }
+
+  /*************** add unmatched stars *************************/
+
+  /** incorporate unmatched image stars, if this star is in field of this tgtcat **/
+  /* these new entries are all written out in UPDATE mode */ 
+  for (off_t i = 0; (i < Nstars) && !options->only_match; i++) {
+    /* make sure there is space for next entry */
+    if (Nmeas >= NMEAS - stars[i].average.Nmeasure) {
+      NMEAS = Nmeas + 1000 + stars[i].average.Nmeasure;
+      REALLOCATE (next_meas, off_t, NMEAS);
+      REALLOCATE (tgtcat[0].measure, Measure, NMEAS);
+    }
+    if (Nave >= NAVE) {
+      NAVE = Nave + 1000;
+      REALLOCATE (tgtcat[0].average, Average, NAVE);
+      if (tgtcat[0].secfilt) {
+	// we only update the secfilt table if it has been allocated for output
+	REALLOCATE (tgtcat[0].secfilt, SecFilt, NAVE*tgtcat[0].Nsecfilt);
+      }
+    }
+
+    if (stars[i].found) continue;
+    if (!IN_REGION (stars[i].average.R, stars[i].average.D)) continue;
+
+    dvo_average_init (&tgtcat[0].average[Nave]);
+    tgtcat[0].average[Nave].R         	   = stars[i].average.R;
+    tgtcat[0].average[Nave].D         	   = stars[i].average.D;
+    tgtcat[0].average[Nave].dR         	   = stars[i].average.dR;
+    tgtcat[0].average[Nave].dD         	   = stars[i].average.dD;
+
+    tgtcat[0].average[Nave].measureOffset  = Nmeas;
+    tgtcat[0].average[Nave].objID     	   = objID;
+    tgtcat[0].average[Nave].catID     	   = catID;
+
+    if (PSPS_ID) {
+      tgtcat[0].average[Nave].extID = CreatePSPSObjectID(tgtcat[0].average[Nave].R, tgtcat[0].average[Nave].D);
+    }
+
+    if (ACCEPT_MOTION) {
+      tgtcat[0].average[Nave].uR    = stars[i].average.uR;
+      tgtcat[0].average[Nave].uD    = stars[i].average.uD;
+      tgtcat[0].average[Nave].duR   = stars[i].average.duR;
+      tgtcat[0].average[Nave].duD   = stars[i].average.duD;
+      tgtcat[0].average[Nave].P     = stars[i].average.P;
+      tgtcat[0].average[Nave].dP    = stars[i].average.dP;
+      tgtcat[0].average[Nave].Tmean = stars[i].average.Tmean;
+    }
+	
+    objID ++;
+
+    // we only update the secfilt table if it has been allocated for output
+    for (int j = 0; tgtcat[0].secfilt && (j < Nsecfilt); j++) {
+      dvo_secfilt_init (&tgtcat[0].secfilt[Nave*Nsecfilt+j], SECFILT_RESET_ALL);
+    }
+
+    // supply the measurements from this detection
+    for (int mSrc = 0; mSrc < stars[i].average.Nmeasure; mSrc++) {
+      dvo_measure_init (&tgtcat[0].measure[Nmeas]);
+      tgtcat[0].measure[Nmeas] = stars[i].measure[mSrc];
+
+      // the following measure elements cannot be set until here:
+      tgtcat[0].measure[Nmeas].dbFlags  = 0;
+      tgtcat[0].measure[Nmeas].averef   = Nave;
+      tgtcat[0].measure[Nmeas].objID    = tgtcat[0].average[Nave].objID;
+      tgtcat[0].measure[Nmeas].catID    = tgtcat[0].catID;
+      
+      /* set the average magnitude if not already set and the photcode.equiv is not 0 */
+      /* in UPDATE mode, this value is not saved; use relphot to recalculate */
+      if (Nsec > -1) { 
+	tgtcat[0].secfilt[Nave*Nsecfilt+Nsec].MpsfChp = PhotCat (&tgtcat[0].measure[Nmeas], MAG_CLASS_PSF);
+      }
+
+      tgtcat[0].average[Nave].Nmeasure ++;
+      Nmeas ++;
+
+      // update the next_meas pointer for this entry (last one for this star is -1)
+      next_meas[Nmeas-1] = (mSrc < stars[i].average.Nmeasure - 1) ? Nmeas : -1;
+    }
+    stars[i].found = TRUE;
+    Nave ++;
+  }
+
+  REALLOCATE (tgtcat[0].average, Average, Nave);
+  REALLOCATE (tgtcat[0].measure, Measure, Nmeas);
+ 
+  if (options->nosort) {
+    tgtcat[0].sorted = FALSE;
+  } else {
+    tgtcat[0].sorted = TRUE;
+    tgtcat[0].measure = sort_measure (tgtcat[0].average, Nave, tgtcat[0].measure, Nmeas, next_meas);
+  }
+
+  /* check if the tgtcat has changed?  if no change, no need to write */
+  tgtcat[0].objID    = objID; // new max value, save on tgtcat close
+  tgtcat[0].Naverage = Nave;
+  tgtcat[0].Nmeasure = Nmeas;
+  tgtcat[0].Nsecfilt_mem = tgtcat[0].secfilt ? Nave*Nsecfilt : 0;
+  if (VERBOSE) fprintf (stderr, "Nstars, Nave, Nmeas: "OFF_T_FMT" "OFF_T_FMT" "OFF_T_FMT", ("OFF_T_FMT" matches)\n",  Nstars,  Nave,  Nmeas, Nmatch);
+
+  free (X1);
+  free (Y1);
+  free (N1);
+  free (X2);
+  free (Y2);
+  free (N2);
+  free (next_meas);
+
+  return (Nmatch);
+}
+
+/* 
+   notes:
+   
+   for finding if a tgtcat star is in an image or an image star is in the tgtcat:
+   
+   tgtcats have boundaries defined by RA and DEC, but they may curve in projection
+   images have boundaries which are lines in pixels coords, but curve in RA and DEC
+   
+   tgtcat[0].found[Ncat] but stars[Nstar].found
+   
+*/
Index: /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs.c
===================================================================
--- /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs.c	(revision 40832)
+++ /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs.c	(revision 40832)
@@ -0,0 +1,44 @@
+# include "addstar.h"
+# include "ukirt_uhs.h"
+
+/* This is the DVO program to upload UKIRT_UHS detections from cvs files supplied by Mike
+   Read into a DVO database.  It is modeled on the loadgaia_dr2 program with column modifications
+   determined by Zhoujian Zhang
+
+   USAGE: loadukirt_uhs -D CATDIR (catdir) (ukirt_uhsfile) [...more files]
+*/
+
+int main (int argc, char **argv) {
+
+  SkyTable *sky;
+  SkyList *skylist = NULL;
+  AddstarClientOptions options;
+
+  // need to construct these options with args_loadtycho...
+  SetSignals ();
+  options = ConfigInit (&argc, argv);
+  options = args_loadukirt_uhs (&argc, argv, options);
+
+  // load the full sky description table:
+  sky = SkyTableLoadOptimal (CATDIR, SKY_TABLE, GSCFILE, TRUE, SKY_DEPTH, VERBOSE);
+  SkyTableSetFilenames (sky, CATDIR, "cpt");
+  
+  // generate the subset matching the user-selected region
+  skylist = SkyListByPatch (sky, -1, &UserPatch);
+
+  int Nstart = 1;
+  int Nend = argc;
+  while (Nstart < Nend) {
+    Nstart = loadukirt_uhs_table (Nstart, Nend, skylist, NULL, argv, &options);
+  }
+
+  SkyTableFree (sky);
+  SkyListFree(skylist);
+  FreePhotcodeTable();
+  free (CATDIR);
+
+  ohana_memcheck (VERBOSE);
+  ohana_memdump (VERBOSE);
+  exit (0);
+}  
+
Index: /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_catalog.c
===================================================================
--- /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_catalog.c	(revision 40832)
+++ /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_catalog.c	(revision 40832)
@@ -0,0 +1,32 @@
+# include "addstar.h"
+# include "ukirt_uhs.h"
+
+int loadukirt_uhs_catalog (UKIRT_UHS_Stars *stars, int Nstars, SkyRegion *region, char *filename, AddstarClientOptions *options) {
+
+  Catalog catalog;
+
+  // now we have all of the loaded stars in this catalog
+  dvo_catalog_init (&catalog, TRUE);
+  catalog.filename = filename;
+  catalog.catformat = dvo_catalog_catformat (CATFORMAT);  // set the default catformat from config data
+  catalog.catmode   = dvo_catalog_catmode (CATMODE);      // set the default catmode from config data
+  catalog.catflags = DVO_LOAD_AVERAGE | DVO_LOAD_SECFILT | DVO_LOAD_MEASURE;
+  catalog.Nsecfilt  = GetPhotcodeNsecfilt ();
+    
+  // an error exit status here is a significant error
+  if (!dvo_catalog_open (&catalog, region, VERBOSE, "w")) {
+    fprintf (stderr, "ERROR: failure to open/create catalog file %s\n", catalog.filename);
+    exit (2);
+  }
+
+  find_matches_ukirt_uhs (region, stars, Nstars, &catalog, options);
+    
+  SetProtect (TRUE);
+  if (!dvo_catalog_save (&catalog, VERBOSE)) { fprintf (stderr, "ERROR: failed to save %s\n", catalog.filename); exit (1); }
+  if (!dvo_catalog_unlock (&catalog)) { fprintf (stderr, "ERROR: failed to unlock %s\n", catalog.filename); exit (1); }
+  SetProtect (FALSE);
+
+  dvo_catalog_free (&catalog);
+
+  return (TRUE);
+}
Index: /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_make_subset.c
===================================================================
--- /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_make_subset.c	(revision 40832)
+++ /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_make_subset.c	(revision 40832)
@@ -0,0 +1,44 @@
+# include "addstar.h"
+# include "ukirt_uhs.h"
+
+// assign stars in the given region to the subset (NOTE: stars are sorted by RA, start is
+// first entry in stars array in this region)
+
+UKIRT_UHS_Stars *loadukirt_uhs_make_subset (UKIRT_UHS_Stars *stars, int Nstars, int start, SkyRegion *region, int *nsubset) {
+
+  int i;
+
+  UKIRT_UHS_Stars *subset = NULL;
+
+  // collect array of (Stars *) stars in a new output catalog
+  int Nsubset = 0;
+  int NSUBSET = 3000;
+  ALLOCATE (subset, UKIRT_UHS_Stars, NSUBSET);
+
+  // find the rest of the stars in this output region
+  for (i = start; i < Nstars; i++) {
+    if (stars[i].flag) continue;
+
+    // check if in skyregion
+    if (stars[i].average.R <  region[0].Rmin) continue;
+    if (stars[i].average.R >= region[0].Rmax) break;
+    if (stars[i].average.D <  region[0].Dmin) continue;
+    if (stars[i].average.D >= region[0].Dmax) continue;
+	  
+    // check if in UserPatch (a GLOBAL)
+    if (stars[i].average.R < UserPatch.Rmin) continue;
+    if (stars[i].average.R > UserPatch.Rmax) break;
+    if (stars[i].average.D < UserPatch.Dmin) continue;
+    if (stars[i].average.D > UserPatch.Dmax) continue;
+	  
+    subset[Nsubset] = stars[i];
+    Nsubset ++;
+
+    stars[i].flag = TRUE;
+
+    CHECK_REALLOCATE (subset, UKIRT_UHS_Stars, NSUBSET, Nsubset, 10000);
+  }
+
+  *nsubset = Nsubset;
+  return subset;
+}
Index: /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_readstars.c
===================================================================
--- /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_readstars.c	(revision 40832)
+++ /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_readstars.c	(revision 40832)
@@ -0,0 +1,369 @@
+# include "addstar.h"
+# include "ukirt_uhs.h"
+
+/* this function reads the values of interest from the UKIRT UHS CSV files:
+
+  data model description: doc/mapping_UKIDSS_unWISE_DVO.pdf
+
+  N -- field (column) number in CSV file (1 counting)
+  | UKIRT field name          : DVO measure field
+  1 sourceID                  : extID
+  2 cuEventID		      : 
+  3 frameSetID		      : imageID (assumes we do not try to go backwards for these photcodes)
+  4 ra			      : R
+  5 dec			      : D
+  6 cx			      : 
+  7 cy			      : 
+  8 cz			      : 
+  9 htmID		      : 
+ 10 l			      : 
+ 11 b			      : 
+ 12 lambda		      : 
+ 13 eta			      : 
+ 14 priOrSec		      : primary bit in dbFlags
+ 15 mergedClassStat	      : 
+ 16 mergedClass		      : psfQFperf or photFlags2
+ 17 pStar		      : psfChisq
+ 18 pGalaxy		      : extNsigma
+ 19 pNoise		      : 
+ 20 pSaturated		      : 
+ 21 eBV			      : 
+ 22 aJ			      : 
+ 23 jHallMag		      : M
+ 24 jHallMagErr		      : dM
+ 25 jPetroMag		      : 
+ 26 jPetroMagErr	      : 
+ 27 jAperMag3		      : Map
+ 28 jAperMag3Err	      : dMap
+ 29 jAperMag4		      : Mkron
+ 30 jAperMag4Err	      : dMron
+ 31 jAperMag6		      : 
+ 32 jAperMag6Err	      : 
+ 33 jGausig		      : 
+ 34 jEll		      : 
+ 35 jPA			      : posangle
+ 36 jErrBits		      : 
+ 37 jDeblend		      : 
+ 38 jClass		      : psfQF
+ 39 jClassStat		      : 
+ 40 jppErrBits		      : photFlags
+ 41 jSeqNum                   : detID
+
+*/
+
+# define iPARSE(NS,NE,FIELD,NAME) {					\
+    cA = iparse_csv_rpt (&ivalue, (NS), (NE), cA, &readStatus);		\
+    if (!readStatus && VERBOSE) {					\
+      gprint (GP_ERR, "suspect field: %d (%s) in %s\n",   (NS), NAME, c0); \
+    }									\
+    lineStatus &= readStatus; FIELD[Nelem] = ivalue; }
+
+# define jPARSE(NS,NE,FIELD,NAME) {					\
+    cA = jparse_csv_rpt (&jvalue, (NS), (NE), cA, &readStatus);		\
+    if (!readStatus && VERBOSE) {					\
+      gprint (GP_ERR, "suspect field: %d (%s) in %s\n",   (NS), NAME, c0); \
+    }									\
+    lineStatus &= readStatus; FIELD[Nelem] = jvalue; }
+
+# define dPARSE(NS,NE,FIELD,NAME) {					\
+    cA = dparse_csv_rpt (&dvalue, (NS), (NE), cA, &readStatus);		\
+    if (!readStatus && VERBOSE) {					\
+      gprint (GP_ERR, "suspect field: %d (%s) in %s\n",   (NS), NAME, c0); \
+    }									\
+    lineStatus &= readStatus; FIELD[Nelem] = dvalue; }
+
+UKIRT_UHS_Stars *loadukirt_uhs_readstars (char *filename, UKIRT_UHS_Stars *stars, int *nstars, AddstarClientOptions *options) {
+
+  int codeJ = GetPhotcodeCodebyName ("UKIRT_UHS_J"); if (!codeJ) Shutdown ("missing photocde UKIRT_UHS_J");
+
+  // XXX I need a UKIRT UHS Mean Epoch
+  time_t UKIRT_UHS_EPOCH = ohana_date_to_sec ("2016/01/01,00:00:00");
+  fprintf (stderr, "WARNING: using an invalid UKIRT_UHS_EPOCH (see loadukirt_uhs_readstars.c:60)\n");
+
+  // read in the full FITS files
+  FILE *f = fopen (filename, "r");
+  if (f == NULL) Shutdown ("can't read ukirt_uhs file: %s", filename);
+
+  int Nskip = 1; // UKIRT UHS CSV files have a single header row (and no special character to mark)
+
+  int Nelem = 0;      // number of valid rows read (vector elements)
+  int NELEM = 10000;  // currently-allocated number of output rows
+
+  // for initial testing, a minimal subset:
+  ALLOCATE_PTR (sourceID, uint64_t, NELEM);
+  ALLOCATE_PTR (Rg,       double, NELEM);
+  ALLOCATE_PTR (Dg,       double, NELEM);
+  ALLOCATE_PTR (chiSq,    double, NELEM);
+  ALLOCATE_PTR (drg,      double, NELEM);
+
+  ALLOCATE_PTR (sourceID, , NELEM);
+  ALLOCATE_PTR (frameSetID,  , NELEM);
+  ALLOCATE_PTR (ra,  , NELEM);
+  ALLOCATE_PTR (dec,  , NELEM);
+  ALLOCATE_PTR (pStar,  , NELEM);
+  ALLOCATE_PTR (pGalaxy,  , NELEM);
+  ALLOCATE_PTR (jHallMag,  , NELEM);
+  ALLOCATE_PTR (jHallMagErr,  , NELEM);
+  ALLOCATE_PTR (jAperMag3,  , NELEM);
+  ALLOCATE_PTR (jAperMag3Err,  , NELEM);
+  ALLOCATE_PTR (jAperMag4,  , NELEM);
+  ALLOCATE_PTR (jAperMag4Err,  , NELEM);
+  ALLOCATE_PTR (jPA,  , NELEM);
+  ALLOCATE_PTR (jClass,  , NELEM);
+  ALLOCATE_PTR (jppErrBits,  , NELEM);
+  ALLOCATE_PTR (photFlags2,  , NELEM);
+  ALLOCATE_PTR (jSeqNum,  , NELEM);
+
+  // we allocate one extra byte into which we never read so there will always be a NULL terminating the string
+  ALLOCATE_PTR (buffer, char, 0x10001);
+  bzero (buffer, 0x10001);
+
+  int Nline_read = 0; // track number of lines read so far (use to skip lines as well)
+
+  // we have a working buffer read from the file. we parse the lines in the working buffer
+  // until we reach the last chunk without an EOL char.  at that point, we shift the start
+  // of the last (partial) line to the start of the buffer and re-fill.
+
+  // we treat \n\r pair as a single EOL char to handle mac files:
+
+  int Nstart = 0; // location of the last valid byte in the buffer (start filling here)
+  int EndOfFile = FALSE;
+  while (!EndOfFile) {
+    int Nbytes = 0x10000 - Nstart;
+    // we have allocated one extra byte into which we never read so there will always be a NULL terminating the string
+    bzero (&buffer[Nstart], Nbytes + 1);
+    int Nread = fread (&buffer[Nstart], 1, Nbytes, f);
+    if (ferror (f)) {
+      perror ("error reading data file");
+      break;
+    }
+    // we still need to parse the rest of the buffer, but there might not be an EOL on the last line
+    if (Nread == 0) {
+      EndOfFile = TRUE;
+    }
+    
+    int bufferStatus = TRUE; 
+    char *c0 = buffer; // c0 always marks the start of a line
+    char *cA = NULL; // cA will carry the curr point within the line
+    while (bufferStatus) {
+
+      // find the end of this current line (\n or \r).  if we hit the end of the buffer (NULL),
+      // attempt to read more data.  finish up when we hit the end of the file
+      char *c1 = strchr (c0, '\n'); // find the end of this current line (also valid for a Mac: \r\n)
+      if (!c1) {
+	c1 = strchr (c0, '\r'); // try \r for Windows files
+      }
+      if (!c1) {
+	Nstart = strlen (c0);
+	if (EndOfFile) {
+	  // if we have reached EOF, we need to do one last pass in case there is a line without a return
+	  c1 = c0 + Nstart;
+	  bufferStatus = FALSE;
+	  if (Nstart == 0) continue; // if we have reached EOF and c0 points at the last valid character, we are done
+	} else {
+	  // if we have not reached EOF, we need to shift the buffer to the start of this line and read more data
+	  memmove (buffer, c0, Nstart);
+	  bufferStatus = FALSE;
+	  continue;
+	}
+      }
+      *c1 = 0; // mark the end of the line 
+      Nline_read ++;
+
+      // skip to the next line (but if EOF, do not overrun buffer)
+      if (Nline_read <= Nskip) { if (!EndOfFile) { c0 = c1 + 1; } continue; }
+
+      // these are not needed: Gaia CSV files do not have any commented-out lines
+      if (*c0 == '#')          { if (!EndOfFile) { c0 = c1 + 1; } continue; }
+      if (*c0 == '!')          { if (!EndOfFile) { c0 = c1 + 1; } continue; }
+
+      // XXX for Gaia, we know which columns we want in advance
+
+      int lineStatus = TRUE;
+      int readStatus;
+      double dvalue;
+      int ivalue;
+      uint64_t jvalue;
+
+      // cA will follow the currently extracted field, c0 points to the start of the line
+      cA = c0;
+
+      // Tref == 2015.5 for DR2
+
+      // the start of the line is the 1st element (fields are 1-counting)
+      jPARSE ( 1,  0, sourceID,     "sourceID");
+      iPARSE ( 3,  1, frameSetID,   "frameSetID");
+      dPARSE ( 4,  3, ra,           "ra");
+      dPARSE ( 5,  4, dec,          "dec");
+      jPARSE (14,  5, priOrSec,     "priOrSec");
+      iPARSE (16, 14, mergedClass,  "mergedClass");
+      dPARSE (17, 16, pStar,        "pStar");
+      dPARSE (18, 17, pGalaxy,      "pGalaxy");
+      dPARSE (23, 18, jHallMag,     "jHallMag");
+      dPARSE (24, 23, jHallMagErr,  "jHallMagErr");
+      dPARSE (27, 24, jAperMag3,    "jAperMag3");
+      dPARSE (28, 27, jAperMag3Err, "jAperMag3Err");
+      dPARSE (29, 28, jAperMag4,    "jAperMag4");
+      dPARSE (30, 29, jAperMag4Err, "jAperMag4Err");
+      dPARSE (35, 30, jPA,          "jPA");
+      dPARSE (38, 35, jClass,       "jClass");
+      iPARSE (40, 38, jppErrBits,   "jppErrBits");
+      iPARSE (41, 40, jSeqNum,      "jSeqNum");
+
+      if (!lineStatus && VERBOSE) {
+	// why do I need to copy temp here, does gprint modify the value of temp?
+	char temp[32];
+	strncpy (temp, c0, 32);
+	temp[31] = 0;
+	gprint (GP_ERR, "skip line %s\n\n", temp);
+      }
+
+      Nelem ++;
+      if (Nelem == NELEM) {
+	NELEM += 1000;
+
+	REALLOCATE (sourceID, uint64_t, NELEM);
+	REALLOCATE (ra, double, NELEM);
+	REALLOCATE (dec, double, NELEM);
+	REALLOCATE (priOrSec, uint64_t, NELEM);
+	REALLOCATE (mergedClass, int, NELEM);
+	REALLOCATE (pStar, double, NELEM);
+	REALLOCATE (pGalaxy, double, NELEM);
+	REALLOCATE (jHallMag, double, NELEM);
+	REALLOCATE (jHallMagErr, double, NELEM);
+	REALLOCATE (jAperMag3, double, NELEM);
+	REALLOCATE (jAperMag3Err, double, NELEM);
+	REALLOCATE (jAperMag4, double, NELEM);
+	REALLOCATE (jAperMag4Err, double, NELEM);
+	REALLOCATE (jPA, double, NELEM);
+	REALLOCATE (jClass,  double, NELEM);
+	REALLOCATE (jppErrBits, int, NELEM);
+	REALLOCATE (jSeqNum, int, NELEM);
+      }
+      if (!EndOfFile) {
+	c0 = c1 + 1;
+      }
+    }
+  }
+
+  // Nelem is now the number of items (objects,stars) read from the Gaia CSV file
+  int NstarsIn = Nelem;
+
+  double Rmin = +360.0;
+  double Rmax = -360.0;
+  double Dmin = +360.0;
+  double Dmax = -360.0;
+
+  // start off where we finished on a previous read
+  int Nstars = *nstars;
+  int NSTARS = Nstars + 0.1*NstarsIn;
+
+  if (!stars) {
+    ALLOCATE (stars, UKIRT_UHS_Stars, NSTARS);
+  } else {
+    REALLOCATE (stars, UKIRT_UHS_Stars, NSTARS);
+  }
+
+  for (int i = 0; i < NstarsIn; i++) {
+
+    Rmin = MIN (Rmin, Rg[i]);
+    Rmax = MAX (Rmax, Rg[i]);
+    Dmin = MIN (Dmin, Dg[i]);
+    Dmax = MAX (Dmax, Dg[i]);
+
+    // each Gaia DR2 source corresponds to 3 measurements: Gg, Gr, Gb:
+
+    dvo_average_init (&stars[Nstars].average);
+    dvo_measure_init (&stars[Nstars].measure);
+
+    stars[Nstars].average.R = ra[i];
+    stars[Nstars].average.D = dec[i];
+    stars[Nstars].average.dR = NAN;
+    stars[Nstars].average.dD = NAN;
+
+    stars[Nstars].flag  = FALSE;
+    stars[Nstars].found = FALSE;
+
+    int isPrimary   = (!priOrSec[i] || (priOrSec[i] == frameSetID[i])) ? ID_MEAS_STACK_PRIMARY : 0x00;
+    int isGalaxy    = (mergedClass[i] == 1)  ? 0x01 : 0x00;
+    int isNoise     = (mergedClass[i] == 0)  ? 0x02 : 0x00;
+    int isStar      = (mergedClass[i] == -1) ? 0x04 : 0x00;
+    int maybeStar   = (mergedClass[i] == -2) ? 0x08 : 0x00;
+    int maybeGalaxy = (mergedClass[i] == -3) ? 0x10 : 0x00;
+    int isSaturated = (mergedClass[i] == -9) ? 0x20 : 0x00;
+
+    int photFlags2 = isGalaxy | isNoise | isStar | maybeStar | maybeGalaxy | isSaturated | isPrimary;
+
+    stars[Nstars].measure[0].extID     = sourceID[i];
+    stars[Nstars].measure[0].imageID   = frameSetID[i];
+    stars[Nstars].measure[0].R         = ra[i];
+    stars[Nstars].measure[0].D         = dec[i];
+    stars[Nstars].measure[0].psfChisq  = pStar[i];
+    stars[Nstars].measure[0].extNsigma = pGalaxy[i];
+    stars[Nstars].measure[0].M         = jHallMag[i];
+    stars[Nstars].measure[0].dM        = jHallMagErr[i];
+    stars[Nstars].measure[0].Map       = jAperMag3[i];
+    stars[Nstars].measure[0].dMap      = jAperMag3Err[i];
+    stars[Nstars].measure[0].Mkron     = jAperMag4[i];
+    stars[Nstars].measure[0].dMkron    = jAperMag4Err[i];
+    stars[Nstars].measure[0].posangle  = jPA[i]; // XXX units
+    stars[Nstars].measure[0].psfQF     = jClass[i];
+
+    stars[Nstars].measure[0].photFlags  = jppErrBits[i];
+    stars[Nstars].measure[0].photFlags2 = photFlags2[i];
+    stars[Nstars].measure[0].detID      = jSeqNum[i];
+
+    stars[Nstars].measure[0].photcode  = codeJ;
+
+    stars[Nstars].average.Nmeasure = 1;
+    Nstars ++;
+
+    CHECK_REALLOCATE (stars, UKIRT_UHS_Stars, NSTARS, Nstars, 10000);
+  }
+
+  free (srcID);
+  free (Rg);
+  free (dRg);
+  free (Dg);
+  free (dDg);
+  free (plx);
+  free (dplx);
+  free (uR);
+  free (duR);
+  free (uD);
+  free (duD);
+  free (nObs);
+  free (nGood);
+  free (goFit);
+  free (chiSq);
+  free (exNoise);
+  free (exNoiseSig);
+  free (nPeriods);
+  free (dupSource);
+  free (gg);
+  free (dgg);
+  free (bg);
+  free (dbg);
+  free (rg);
+  free (drg);
+  free (procMode);
+
+  free (buffer);
+
+  *nstars = Nstars;
+  return (stars);
+}
+
+int loadukirt_uhs_sortStars (UKIRT_UHS_Stars *stars, int Nstars) {
+
+# define SWAPFUNC(A,B){ UKIRT_UHS_Stars temp = stars[A]; stars[A] = stars[B]; stars[B] = temp; }
+# define COMPARE(A,B)(stars[A].average.R < stars[B].average.R)
+
+  OHANA_SORT (Nstars, COMPARE, SWAPFUNC);
+
+# undef SWAPFUNC
+# undef COMPARE
+  
+  return TRUE;
+}
+
Index: /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_table.c
===================================================================
--- /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_table.c	(revision 40832)
+++ /branches/eam_branches/ohana.20190329/src/addstar/src/loadukirt_uhs_table.c	(revision 40832)
@@ -0,0 +1,64 @@
+# include "addstar.h"
+# include "ukirt_uhs.h"
+# define NSTARS_MAX 1000000
+
+int loadukirt_uhs_table (int Nstart, int Nend, SkyList *skylistInput, HostTable *hosts, char **filename, AddstarClientOptions *options) {
+  OHANA_UNUSED_PARAM(hosts);
+  
+  int i;
+
+  int Nstars = 0;
+  UKIRT_UHS_Stars *stars = NULL;
+  for (i = Nstart; (Nstars < NSTARS_MAX) && (i < Nend); i++) {
+    // read the next file and append to the current set of stars
+    fprintf (stderr, "loading %s\n", filename[i]);
+    stars = loadukirt_uhs_readstars (filename[i], stars, &Nstars, options);
+  }
+  Nstart = i; // we pass back the entry for the next file to be read
+  if (!stars) return Nstart;
+
+  fprintf (stderr, "writing %d stars to dvo\n", Nstars);
+
+  // sort the stars by RA
+  loadukirt_uhs_sortStars (stars, Nstars);
+
+  // scan through the stars, loading the containing catalogs
+  // skip through table for unsaved stars
+  for (i = 0; i < Nstars; i++) {
+    if (stars[i].flag) continue;
+
+    // scan forward until we read the UserPatch
+    if (stars[i].average.R < UserPatch.Rmin) continue;
+    if (stars[i].average.R > UserPatch.Rmax) break;
+    if (stars[i].average.D < UserPatch.Dmin) continue;
+    if (stars[i].average.D > UserPatch.Dmax) continue;
+
+    // identify the relevant catalog
+    SkyList *skylist = SkyRegionByPoint_List (skylistInput, -1, stars[i].average.R, stars[i].average.D);
+    if (skylist[0].Nregions == 0) {
+      SkyListFree (skylist);
+      continue;
+    }
+    SkyRegion *region = skylist[0].regions[0];
+
+    // select stars matching this region
+    int Nsubset;
+    UKIRT_UHS_Stars *subset = loadukirt_uhs_make_subset (stars, Nstars, i, region, &Nsubset);
+
+    // In parallel mode, write out the subset to a disk file.  Block until a remote host
+    // is available.  In serial mode, just match against the appropriate region and save
+    // NOTE: disable parallel mode for now: 
+    // loadukirt_uhs_save_remote (subset, Nsubset, hosts, region, skylist[0].filename[0], options);
+    loadukirt_uhs_catalog (subset, Nsubset, region, skylist[0].filename[0], options);
+    free (subset);
+    SkyListFree (skylist);
+  }
+
+  // wait for last remote clients to finish
+  // NOTE: disable parallel mode for now: 
+  // harvest_all ();
+
+  free (stars);
+  return Nstart;
+}
+
