Changeset 36466

Timestamp:

Feb 1, 2014, 6:04:48 AM (12 years ago)

Author:

eugene

Message:

extending boundary tree stuff to support local tessallations (not ready yet)

Location:

branches/eam_branches/ipp-20131211/Ohana/src/relphot/src

Files:

• BoundaryTreeOps.c (modified) (1 diff)
• ImageOps.c (modified) (1 diff)
• args.c (modified) (1 diff)
• setMrelCatalog.c (modified) (1 diff)

branches/eam_branches/ipp-20131211/Ohana/src/relphot/src/BoundaryTreeOps.c

@@ -1 +1 @@
 # include "relphot.h"
 
-// XXX for the moment, only load one boundary tree at a time
-// XXX in fact, only allow RINGS.V3...
+static int Nfail = 0;
+static int Ntess = 0;
+static TessllationBoundary *tess = NULL;
 
-static BoundaryTree *tree = NULL;
-
-int BoundaryTreePrimaryCell (char *primaryCellName, double ra, double dec) {
+// for the given ra,dec : find the valid tessallation and the containing projection / skycell IDs
+int TessallationPrimaryCellIDs (int *tessID, int *projID, int *skycellID, double ra, double dec) {
 
   int zone, band;
 
-  if (!primaryCellName) return FALSE;
-
-  primaryCellName[0] = 0;
-
-  if (!tree) return FALSE;
-
-  if (!BoundaryTreeCellCoords (tree, &zone, &band, ra, dec)) {
-    fprintf (stderr, "mismatch!\n");
-    return FALSE;
-  }
-  
-# define USE_PROJECTION_CELL 0
-# if (USE_PROJECTION_CELL)
-  snprintf (primaryCellName, DVO_MAX_PATH, "RINGS.V3.%s", tree->name[zone][band]);
-# else
-
-  // I have ra, dec, and the primary projection cell.  In order to choose the primary skycell,
-  // I just need to project to ra,dec to X,Y based on the center of the cell and then get the subdivision right.
-    
-  double x = 0.0;
-  double y = 0.0;
-  BoundaryTreeProjection (&x, &y, ra, dec, tree, zone, band);
-  
-  int xi = x / tree->dX[zone][band];
-  int yi = y / tree->dY[zone][band];
-  int N = xi + tree->NX_SUB * yi;
-  
-  // XXX short-circuit this for now (we should use this code if we make tree have more variable NX,NY values
-  // char format[24], skycellname[128];
-  // int Ndigit = (int)(log10(tree->NX_SUB*tree->NY_SUB)) + 1 ;
-  // snprintf (format, 24, "%s.%%0%dd", tree->name[zone][band], Ndigit);
-  // snprintf (skycellname, 128, format, N);
-
-  snprintf (primaryCellName, DVO_MAX_PATH, "RINGS.V3.%s.%03d", tree->name[zone][band], N);
-# endif
-
-  return TRUE;
-}
-
-int BoundaryTreePrimaryCellIDs (int *projID, int *skycellID, double ra, double dec) {
-
-  int zone, band;
-
+  *tessID = -1;
   *projID = -1;
   *skycellID = -1;
 
-  if (!tree) return FALSE;
+  if (!tess) return FALSE;
 
-  if (!BoundaryTreeCellCoords (tree, &zone, &band, ra, dec)) {
-    fprintf (stderr, "mismatch!\n");
+  // find the tessellation which includes this location
+  int myTess = -1;
+  for (i = 0; i < Ntess; i++) {
+    if (ra  <  tess[i].Rmin) continue;
+    if (ra  >= tess[i].Rmax) continue;
+    if (dec <  tess[i].Dmin) continue;
+    if (dec >= tess[i].Dmax) continue;
+    myTess = i;
+    break;
+  }
+  if (myTess < 0) {
+    if (Nfail < 100) {
+      fprintf (stderr, "no matching tessallation @ %f,%f\n", ra, dec);
+    }
+    Nfail ++;
     return FALSE;
   }
+
+  if (tess[myTess].type == TESS_NONE) {
+    if (Nfail < 100) {
+      fprintf (stderr, "tess has invalid type\n");
+      Nfail ++;
+    }
+    return FALSE;
+  }
+
+  if (tess[myTess].type == TESS_LOCAL) {
+    // I have ra, dec, and the primary projection cell.  In order to choose the primary skycell,
+    // I just need to project to ra,dec to X,Y based on the center of the cell and then get the subdivision right.
+    
+    double x = 0.0;
+    double y = 0.0;
+    TessallationLocalProjection (&x, &y, ra, dec, &tess[myTess]);
   
-  // I have ra, dec, and the primary projection cell.  In order to choose the primary skycell,
-  // I just need to project to ra,dec to X,Y based on the center of the cell and then get the subdivision right.
+    int xi = x / tess[myTess].dX;
+    int yi = y / tess[myTess].dY;
+    int N = xi + tess[myTess].NX_SUB * yi;
+  
+    *tessID = myTess;
+    *projID = tess[myTess].projID;
+    *skycellID = N;
+
+    return TRUE;
+  }
+
+  if (tess[myTess].type == TESS_RINGS) {
+    if (!BoundaryTreeCellCoords (tess[myTess].tree, &zone, &band, ra, dec)) {
+      fprintf (stderr, "mismatch!\n");
+      return FALSE;
+    }
+  
+    // I have ra, dec, and the primary projection cell.  In order to choose the primary skycell,
+    // I just need to project to ra,dec to X,Y based on the center of the cell and then get the subdivision right.
     
-  double x = 0.0;
-  double y = 0.0;
-  BoundaryTreeProjection (&x, &y, ra, dec, tree, zone, band);
+    double x = 0.0;
+    double y = 0.0;
+    BoundaryTreeProjection (&x, &y, ra, dec, tess[myTess].tree, zone, band);
   
-  int xi = x / tree->dX[zone][band];
-  int yi = y / tree->dY[zone][band];
-  int N = xi + tree->NX_SUB * yi;
+    int xi = x / tess[myTess].tree->dX[zone][band];
+    int yi = y / tess[myTess].tree->dY[zone][band];
+    int N = xi + tess[myTess].tree->NX_SUB * yi;
   
-  *projID = tree->projID[zone][band];
-  *skycellID = N;
+    *tessID = myTess;
+    *projID = tess[myTess].tree->projID[zone][band];
+    *skycellID = N;
 
-  return TRUE;
+    return TRUE;
+  }
 }
 
-int load_tree (char *treefile) {
+int TessallationIDsByImageName (int *tessID, int *projID, int *skycellID, char *name) {
 
-  tree = BoundaryTreeLoad (treefile);
-  if (!tree) {
-    fprintf (stderr, "failed to load boundary tree %s\n", treefile);
+  *tessID = -1;
+  *projID = -1;
+  *skycellID = -1;
+
+  if (!tess) return FALSE;
+
+  for (i = 0; i < Ntess; i++) {
+    // do this with a strhash of some kind?
+    if (!strncmp (name, tess[i].basename, tess[i].Nbasename)) {
+      *tessID    = i;
+      *projID    = atoi(&name[tess[i].projectIDoff]);
+      *skycellID = atoi(&name[tess[i].skycellIDoff]);
+      return TRUE;
+    }
+  }
+  return FALSE;
+}
+
+int load_tess (char *tessfile) {
+
+  tess = TessBoundaryLoad (tessfile, &Ntess);
+  if (!tess) {
+    fprintf (stderr, "failed to load tessallation boundary file %s\n", tessfile);
     exit (2);
   }

branches/eam_branches/ipp-20131211/Ohana/src/relphot/src/ImageOps.c

@@ -87 +87 @@
     imageIdx[i] = i;
     imageIDs[i] = image[i].imageID;
+
+    tessID[i] = -1;
     projectID[i] = -1;
     skycellID[i] = -1;
-    if (!strncmp (image[i].name, "RINGS.V3.skycell", strlen("RINGS.V3.skycell"))) {
-      projectID[i] = atoi(&image[i].name[17]);
-      skycellID[i] = atoi(&image[i].name[22]);
-    }
+
+    TessallationIDsByImageName (&tessID[i], &projectID[i], &skycellID[i], image[i].name);
   }
 
   // sort the image index by the IDs
+  // XXX does this break the imageID <-> projectID, etc match?
   llsortpair (imageIDs, imageIdx, Nimage);
 }

branches/eam_branches/ipp-20131211/Ohana/src/relphot/src/args.c

@@ -217 +217 @@
   }
 
+  // XXX should we load a tree from CATDIR by default?
+  // NOTE: a given catdir needs an appropriate boundary tree
   BOUNDARY_TREE = NULL;
   if ((N = get_argument (argc, argv, "-boundary-tree"))) {

branches/eam_branches/ipp-20131211/Ohana/src/relphot/src/setMrelCatalog.c

@@ -405 +405 @@
   float Msys = 0, Mcal= 0, Mmos = 0, Mgrid = 0;
 
-  // set the name of the primary skycell (this is used in a strcmp to match the skycells in stack detections)
-  // XXX drop BoundaryTreePrimaryCell(primaryCell, average[0].R, average[0].D);
+  // set the primary projection cell and skycell for this coordinate
+  int tessID    = -1;
   int projectID = -1;
   int skycellID = -1;
-  BoundaryTreePrimaryCellIDs(&projectID, &skycellID, average[0].R, average[0].D);
+  TessallationTreePrimaryCellIDs(&tessID, &projectID, &skycellID, average[0].R, average[0].D);
 
   int NstackGood = 0;