Changeset 40528

Timestamp:

Sep 25, 2018, 10:21:11 AM (8 years ago)

Author:

eugene

Message:

added CFIS tessellation to skycells

Location:

trunk/Ohana/src/addstar

Files:

• include/skycells.h (modified) (3 diffs)
• src/args_skycells.c (modified) (5 diffs)
• src/sky_tessalation.c (modified) (3 diffs)

trunk/Ohana/src/addstar/include/skycells.h

@@ -12 +12 @@
 # include <glob.h>
 
-enum {SQUARES, TRIANGLES, LOCAL, RINGS, TAMAS};
+enum {SQUARES, TRIANGLES, LOCAL, RINGS, TAMAS, CFIS};
 enum {TETRAHEDRON, CUBE, OCTOHEDRON, DODECAHEDRON, ICOSAHEDRON};
 
@@ -96 +96 @@
 int          sky_tessellation_rings         PROTO((FITS_DB *db, int level, int Nmax));
 int          sky_tessellation_tamas         PROTO((FITS_DB *db, int level, int Nmax));
+int          sky_tessellation_cfis          PROTO((FITS_DB *db, int level, int Nmax));
 
 int          sky_triangle_to_image          PROTO((Image *image, SkyTriangle *triangle));
@@ -106 +107 @@
 SkyRectangle *sky_rectangle_ring            PROTO((float dec, float dDEC, int *nring, char *format));
 SkyRectangle *sky_rectangle_tamas           PROTO((double *Dec, double dm, double halfa, double halftheta, int *nring, char *format));
+SkyRectangle *sky_rectangle_cfis            PROTO((double dec, int *nring, char *format));
 
 SkyTriangle *sky_divide_triangles           PROTO((SkyTriangle *in, int *ntriangles));

trunk/Ohana/src/addstar/src/args_skycells.c

@@ -40 +40 @@
     if (!strcasecmp (argv[N], "tamas")) {
       MODE = TAMAS;
+    }
+    if (!strcasecmp (argv[N], "cfis")) {
+      MODE = CFIS;
     }
     remove_argument (N, &argc, argv);
@@ -115 +118 @@
   /* pixel scale (arcsec/pixel) */
   SCALE = 1.0;
+  if (MODE == CFIS) {
+    SCALE = 0.185768447409; // default CFIS value
+  }
   if ((N = get_argument (argc, argv, "-scale"))) {
     remove_argument (N, &argc, argv);
@@ -174 +180 @@
 
   CELLSIZE = 4.0;
-  if ((MODE == RINGS) && (N = get_argument (argc, argv, "-cellsize"))) {
+  if (MODE == TAMAS) CELLSIZE = 3.955;
+  if (MODE == CFIS) CELLSIZE = 0.5;
+
+  int checkCellsize = (MODE == RINGS) || (MODE == TAMAS) || (MODE == CFIS);
+
+  if (checkCellsize && (N = get_argument (argc, argv, "-cellsize"))) {
     remove_argument (N, &argc, argv);
     CELLSIZE = strtod (argv[N], &ptr);
@@ -183 +194 @@
     remove_argument (N, &argc, argv);
   }
-  if (MODE == TAMAS) {
-    CELLSIZE = 3.955;
-    if ((N = get_argument (argc, argv, "-cellsize"))) {
-      remove_argument (N, &argc, argv);
-      CELLSIZE = strtod (argv[N], &ptr);
-      if ((*ptr != 0) || (CELLSIZE < 0.0)) {
-        fprintf (stderr, "-cellsize requires a floating-point argument\n");
-        help ();
-      }  
-      remove_argument (N, &argc, argv);
-    }
-  }
 
   OVERLAP_RA = 0;
@@ -230 +229 @@
   fprintf (stderr, "     TRIANGLE                 : generate triangular skycells using icosahedron base solid\n");
   fprintf (stderr, "     RINGS                    : generate rectangular skycells using declination strips\n");
+  fprintf (stderr, "     CFIS                     : generate rectangular skycells using CFIS declination strips\n");
   fprintf (stderr, "     LOCAL                    : generate a local tessellation around a spot on the sky\n");
   fprintf (stderr, "                                 (Note that this tessellation does not cover the full sky)\n");

trunk/Ohana/src/addstar/src/sky_tessalation.c

@@ -26 +26 @@
       sky_tessellation_tamas (db, level, Nmax);
       return TRUE;
+    case CFIS:
+      sky_tessellation_cfis (db, level, Nmax);
+      return TRUE;
     default:
       break;
@@ -235 +238 @@
 
   free (image);
+  return (TRUE);
+}
+
+// the CADC / CFIS tessellation uses the constant DEC offsets and simple cos(DEC)-scaled 
+// RA offsets.  The skycell names are defined to match the RA,DEC sequence
+int sky_tessellation_cfis (FITS_DB *db, int level, int Nmax) {
+  OHANA_UNUSED_PARAM(level);
+  OHANA_UNUSED_PARAM(Nmax);
+
+  char format[32];
+
+  int Ndec = 0;
+
+  // generate the a collection of rectangles for each ring
+  for (double dec = -90.0; dec < +90.01; dec += CELLSIZE, Ndec ++) {
+
+    snprintf (format, 32, "skycell.%%03d.%03d", Ndec);
+
+    int Nring;
+    SkyRectangle *ring = sky_rectangle_cfis (dec, &Nring, format);
+    if (!ring) continue;
+
+    // nominal CFIS uses a single skycell per projection center
+
+    // subdivide each image (Nx x Ny subcells)
+    int Nimage = NX_SUB*NY_SUB*Nring;
+    ALLOCATE_PTR (image, Image, Nimage);
+    for (int j = 0; j < Nring; j++) {
+      // convert the SkyRectangles to Images for output
+      sky_subdivide_image (&image[j*NX_SUB*NY_SUB], &ring[j], NX_SUB, NY_SUB);
+      // printf("%s %8.2f %8.2f\n", ring[j].name, ring[j].coords.crval1, ring[j].coords.crval2);
+    }
+
+    /* add the new images and save */
+    dvo_image_addrows (db, image, Nimage);
+    SetProtect (TRUE);
+    dvo_image_update (db, VERBOSE);
+    SetProtect (FALSE);
+    dvo_image_clear_vtable (db);
+    
+    free (ring);
+    free (image);
+  }    
   return (TRUE);
 }
@@ -623 +669 @@
 
   return (TRUE);
+}
+
+// define the parameters of a projection centers for this DEC band 
+// dec : ~ center of band in Dec
+// dDEC : approximate height
+// nring : number of cells generated for this ring
+// format : guide to generate the filenames (c-type string format)
+
+// for skycell xxx.yyy, Dcenter = (yyy / 2 - 900), Rcenter = (xxx / 2 / dcos(Dcenter))
+
+SkyRectangle *sky_rectangle_cfis (double dec, int *nring, char *format) {
+
+  int i;
+
+  // Subdivide the DEC center line into an integer number of segments:
+  float dRA = CELLSIZE / cos(RAD_DEG*dec); // dRA is a size in RA degrees == \alpha_n
+  int nRA = (360.0 + 2*CELLSIZE) / dRA;      // CELLSIZE is projection center spacing (DEC-direction)
+
+  ALLOCATE_PTR (ring, SkyRectangle, nRA);
+
+  for (i = 0; i < nRA; i++) {
+    memset (&ring[i], 0, sizeof(SkyRectangle));
+
+    InitCoords (&ring[i].coords, "DEC--TAN");
+    ring[i].coords.crval1 = i*dRA;
+    ring[i].coords.crval2 = dec;
+
+    ring[i].coords.pc1_1 = +1.0 * X_PARITY;
+  
+    // we need to make this a user-defined option
+    // for the default CFIS cells, PADDING is 0 (hard-wired into NX,NY = 10000)
+    // ring[i].NX = NX*(1.0 + PADDING);
+    // ring[i].NY = NY*(1.0 + PADDING);
+    ring[i].NX = 10000;
+    ring[i].NY = 10000;
+    ring[i].photcode = 1; // this needs to be set more sensibly
+
+    // crpix1,crpix2 is the projection center
+    ring[i].coords.crpix1 = 0.5*ring[i].NX;
+    ring[i].coords.crpix2 = 0.5*ring[i].NY;
+
+    // user-supplied pixel scale
+    ring[i].coords.cdelt1 = SCALE / 3600.0;
+    ring[i].coords.cdelt2 = SCALE / 3600.0;
+
+    // CFIS uses names with centers based on ra,dec
+    // format = skycell.%03d.%03d
+    snprintf (ring[i].name, DVO_IMAGE_NAME_LEN, format, i);
+  }
+
+  *nring = nRA;
+  return (ring);
 }