Changeset 27594 for trunk

Timestamp:

Apr 4, 2010, 9:54:36 AM (16 years ago)

Author:

eugene

Message:

adding function avmatch to extract average fields based on a coordinate match

Location:

trunk/Ohana/src/opihi

Files:

• dvo/Makefile (modified) (2 diffs)
• dvo/avextract.c (modified) (1 diff)
• dvo/avmatch.c (added)
• dvo/dbCmdlineFields.c (modified) (1 diff)
• dvo/find_matches.c (added)
• dvo/init.c (modified) (2 diffs)
• dvo/mextract.c (modified) (1 diff)
• dvo/skyregion.c (modified) (1 diff)
• include/dvoshell.h (modified) (1 diff)

trunk/Ohana/src/opihi/dvo/Makefile

@@ -28 +28 @@
 $(SRC)/dvomisc.$(ARCH).o                \
 $(SRC)/region_list.$(ARCH).o            \
+$(SRC)/find_matches.$(ARCH).o           \
 $(SRC)/photometry.$(ARCH).o             \
 $(SRC)/dbBooleanCond.$(ARCH).o          \
@@ -43 +44 @@
 cmds = \
 $(SRC)/avextract.$(ARCH).o              \
+$(SRC)/avmatch.$(ARCH).o                \
 $(SRC)/badimages.$(ARCH).o              \
 $(SRC)/calextract.$(ARCH).o             \

trunk/Ohana/src/opihi/dvo/avextract.c

@@ -182 +182 @@
 
   if ((argc > N + 1) && !strcasecmp (argv[N+1], "fields")) {
-    gprint (GP_ERR, " USAGE: avextract field[,field,field...] where (expression)\n");
     gprint (GP_ERR, "  RA : right ascension (J2000) [degrees]\n");
     gprint (GP_ERR, "  DEC : declination [degrees]\n");

trunk/Ohana/src/opihi/dvo/dbCmdlineFields.c

@@ -57 +57 @@
   status = FALSE;
 
-  // examine each argv[i] entry until we reach a 'where' or a 'matched' 
+  // examine each argv[i] entry until we reach a 'where', a 'matched', or the end of the line
   for (i = 1; (i < argc) && strcasecmp (argv[i], "where") && strcasecmp (argv[i], "match"); i++) {
     // split the word by ","

trunk/Ohana/src/opihi/dvo/init.c

@@ -2 +2 @@
 
 int avextract       PROTO((int, char **));
+int avmatch         PROTO((int, char **));
 int badimages       PROTO((int, char **));
 int calextract      PROTO((int, char **));
@@ -55 +56 @@
 static Command cmds[] = {  
   {1, "avextract",   avextract,    "extract average data values"},
+  {1, "avmatch",     avmatch,      "extract average data values matched to RA,DEC points"},
   {1, "badimages",   badimages,    "look for images with anomalous astrometry"},
   {1, "calextract",  calextract,   "extract photometry calibration"},

trunk/Ohana/src/opihi/dvo/mextract.c

@@ -202 +202 @@
 
   if ((argc > N + 1) && !strcasecmp (argv[N+1], "fields")) {
-    gprint (GP_ERR, " USAGE: avextract field[,field,field...] where (expression)\n");
     gprint (GP_ERR, "  RA : right ascension (J2000) for detection [degrees]\n");
     gprint (GP_ERR, "  DEC : declination for detection [degrees]\n");

trunk/Ohana/src/opihi/dvo/skyregion.c

@@ -44 +44 @@
   return TRUE;
 }
+
+/* find region which overlaps c at given depth (-1 : populated ) */
+int SkyRegionByPoint_r (SkyTable *table, SkyList *list, int depth, double ra, double dec) {
+  
+  int i, Ns, Ne, No;
+  SkyRegion *skyregion;
+
+  list[0].Nregions = 0;
+  list[0].ownElements = FALSE; // this list is only holding a view to the elements
+
+  skyregion = table[0].regions;
+
+  Ns = 0;
+  Ne = 1;
+
+  while (1) {
+    No = -1;
+    for (i = Ns; (i < Ne) && (i < table[0].Nregions); i++) {
+      if (ra  < skyregion[i].Rmin) continue;
+      if (ra  > skyregion[i].Rmax) continue;
+      if (dec < skyregion[i].Dmin) continue;
+      if (dec > skyregion[i].Dmax) continue;
+      No = i;
+      break;
+    }
+    if (No == -1) return (FALSE);
+    if ((depth == -1) && (skyregion[No].table)) break;
+    if (depth == skyregion[No].depth) break;
+    if ((depth > skyregion[No].depth) && !skyregion[No].child) return (FALSE);
+
+    /* need to check Ns, Ne, or guarantee valid range */
+    Ns = skyregion[No].childS;
+    Ne = skyregion[No].childE;
+  }
+
+  list[0].regions[0] = &skyregion[No];
+  list[0].filename[0] = table[0].filename[No];
+  list[0].Nregions = 1;
+  return (TRUE);
+}
+
+SkyList *SelectRegionsByCoordVectors (Vector *RA, Vector *DEC) {
+  
+  int i, j, Npts, Nout, NOUT, *found;
+  double ra, dec;
+  SkyList *new, *list;
+  SkyTable *sky;
+
+  sky = GetSkyTable ();
+
+  Npts = RA->Nelements;
+  ALLOCATE (found, int, Npts);
+  memset (found, 0, Npts*sizeof(int));
+
+  ALLOCATE (new, SkyList, 1);
+  ALLOCATE (new[0].regions,  SkyRegion *, 1);
+  ALLOCATE (new[0].filename,  char *, 1);
+
+  // output list
+  Nout = 0;
+  NOUT = 100;
+  ALLOCATE (list, SkyList, 1);
+  ALLOCATE (list[0].regions,  SkyRegion *, NOUT);
+  ALLOCATE (list[0].filename,  char *, NOUT);
+  list[0].Nregions = 0;
+  list[0].ownElements = FALSE; // this list is only holding a view to the elements
+
+  for (i = 0; i < Npts; i++) {
+    if (found[i]) continue;
+
+    SkyRegionByPoint_r (sky, new, -1, RA->elements.Flt[i], DEC->elements.Flt[i]);
+    if (new->Nregions == 0) continue;
+    assert (new->Nregions == 1);
+    // append new region to output list
+    list[0].regions[Nout] = new[0].regions[0];
+    list[0].filename[Nout] = new[0].filename[0];
+    Nout ++;
+    if (Nout >= NOUT) {
+        NOUT += 100;
+        REALLOCATE (list[0].regions, SkyRegion *, NOUT);
+        REALLOCATE (list[0].filename, char *, NOUT);
+    }
+    found[i] = TRUE;
+
+    // scan over the remaining points to find any that lie within this region
+    // if we sorted the ra,dec inputs we could break out of this more quickly
+    for (j = i + 1; j < Npts; j++) {
+      ra = RA->elements.Flt[j];
+      dec = DEC->elements.Flt[j];
+      if (ra  < new[0].regions[0][0].Rmin) continue;
+      if (ra  > new[0].regions[0][0].Rmax) continue;
+      if (dec < new[0].regions[0][0].Dmin) continue;
+      if (dec > new[0].regions[0][0].Dmax) continue;
+      found[j] = TRUE;
+    }
+  }
+
+  free (found);
+  SkyListFree (new);
+
+  return (list);
+}
+

trunk/Ohana/src/opihi/include/dvoshell.h

@@ -314 +314 @@
 SkyRegionSelection *SetRegionSelection    PROTO((int *argc, char **argv));
 
+int           SkyRegionByPoint_r    PROTO((SkyTable *table, SkyList *list, int depth, double ra, double dec));
+SkyList      *SelectRegionsByCoordVectors PROTO((Vector *RA, Vector *DEC));
+
+int           find_matches_by_vectors PROTO((SkyRegion *region, Catalog *catalog, Vector *RAvec, Vector *DECvec, float RADIUS, off_t *index));
+
 int           SetImageSelection     PROTO((int mode, SkyRegionSelection *selection));
 int           SetPhotSelections     PROTO((int *argc, char **argv, int Nparams));