Changeset 30783

Timestamp:

Mar 3, 2011, 2:38:35 PM (15 years ago)

Author:

bills

Message:

integrate fox for the problem reported in ticket 1381

File:

• tags/ipp-20110218/Ohana/src/getstar/src/MatchCoords.c (modified) (6 diffs, 1 prop)

tags/ipp-20110218/Ohana/src/getstar/src/MatchCoords.c

@@ -3 +3 @@
 # define FLT_MAX 1e32
 # endif
+int ComputeLMLimits(double xmin, double xmax, double ymin, double ymax, 
+            double *Lmin, double *Lmax, double *Mmin, double *Mmax, Coords *coords);
 
 /* given coordinate, find images in list that contain the point */
@@ -12 +14 @@
   double Xi[4], Yi[4];  /* image and original corners */
   double xmin, xmax, ymin, ymax;
+  double Lmin, Lmax, Mmin, Mmax;
 
 
@@ -20 +23 @@
      also define the vtable entries for the images we keep  */
 
+  int haveLimits = 0;
   for (i = 0; i < NdbImages; i++) {
+      haveLimits = 0;
 
     if (!WITH_PHU && !strcmp (&dbImages[i].coords.ctype[4], "-DIS")) continue;
@@ -30 +35 @@
     /* define image corners */
     SetImageCorners (Xi, Yi, &dbImages[i]);
-    // Xi[4] = Xi[0]; Yi[4] = Yi[0];
+    // Xi[3] = Xi[0]; Yi[3] = Yi[0];
 
     ymin = xmin = +FLT_MAX;
@@ -52 +57 @@
             continue;
         }
-
         
-        if (VERBOSE) fprintf(stderr, "%ld %8.1f %8.1f %s\n", i, x, y, dbImages[i].name);
+       if (VERBOSE) fprintf(stderr, "%ld %8.1f %8.1f %s\n", i, x, y, dbImages[i].name);
 
         if ((x >= xmin && x <= xmax) && (y >= ymin && y <= ymax)) {
+            if (VERBOSE) fprintf(stderr, "\tpotential overlap for point %d image %i\n", j, (int) i);
+
+            // We may be using the astrometry far outside it's range of validity.
+            // As a sanity check make sure that the focal plane coords are in bounds.
+            // See ticket 1381
+            if (!haveLimits) {
+                // We defer computing the focal plane limits until we need them
+                status = ComputeLMLimits(xmin, xmax, ymin, ymax, 
+                    &Lmin, &Lmax, &Mmin, &Mmax, &dbImages[i].coords);
+                if (!status) {
+                    continue;
+                }
+                haveLimits = 1;
+            }
+
+            double L, M;
+            status = RD_to_LM(&L, &M, pt->ra, pt->dec, &dbImages[i].coords);
+            if (!status) {
+                // This won't happen because RD_to_XY succeeded above
+                continue;
+            }
+
+            if (VERBOSE) fprintf(stderr, "L: %f M: %f\n", L, M);
+
+            if (L < Lmin || L > Lmax || M < Mmin || M > Mmax) {
+                if (VERBOSE) fprintf(stderr, "\trejecting match for because focal plane coordinates are out of range\n");
+                continue;
+            }
+            if (VERBOSE) fprintf(stderr, "\tmatch confirmed\n");
+
             totalMatches++;
 
@@ -88 +122 @@
   }
 }
+
+// compute bounding box of this image in focal plane coordinates
+int ComputeLMLimits(double xmin, double xmax, double ymin, double ymax, 
+            double *Lmin, double *Lmax, double *Mmin, double *Mmax, Coords *coords)
+{
+    double L0, L1, M0, M1;
+    int status = XY_to_LM(&L0, &M0, xmin, ymin, coords);
+    if (!status) {
+        fprintf(stderr, "failed to transform image corner to LM\n");
+        return 0;
+    }
+    status = XY_to_LM(&L1, &M1, xmax, ymax, coords);
+    if (!status) {
+        fprintf(stderr, "failed to transform image corner to LM\n");
+        return 0;
+    }
+
+    if (L1 > L0) {
+        *Lmin = L0;
+        *Lmax = L1;
+    } else {
+        *Lmin = L1;
+        *Lmax = L0;
+    }
+    if (M1 > M0) {
+        *Mmin = M0;
+        *Mmax = M1;
+    } else {
+        *Mmin = M1;
+        *Mmax = M0;
+    }
+    if (VERBOSE) fprintf(stderr, "Lmin: %f Lmax: %f Mmin: %f Mmax: %f\n", *Lmin, *Lmax, *Mmin, *Mmax);
+    return 1;
+}