Index: /trunk/Ohana/src/addstar/include/skycells.h =================================================================== --- /trunk/Ohana/src/addstar/include/skycells.h (revision 12772) +++ /trunk/Ohana/src/addstar/include/skycells.h (revision 12773) @@ -11,4 +11,6 @@ # include # include + +enum {SQUARES, TRIANGLES}; typedef struct { @@ -49,4 +51,5 @@ int VERBOSE; +int MODE; void SetProtect PROTO((int mode)); @@ -57,5 +60,4 @@ int ConfigInit_skycells (int *argc, char **argv); -Image *sky_tessalation (int level, int *nimages); int sky_triangle_to_image (Image *image, SkyTriangle *triangle); int sky_triangle_coords (SkyTriangle *triangle); @@ -65,2 +67,6 @@ int sky_tessalation_init (); int sky_triangle_to_rectangle (Image *image, SkyTriangle *triangle); + +Image *sky_tessalation (int level, int *nimages, int mode); +Image *sky_tessalation_triangles (int level, int *nimages); +Image *sky_tessalation_squares (int level, int *nimages); Index: /trunk/Ohana/src/addstar/src/args_skycells.c =================================================================== --- /trunk/Ohana/src/addstar/src/args_skycells.c (revision 12772) +++ /trunk/Ohana/src/addstar/src/args_skycells.c (revision 12773) @@ -19,7 +19,14 @@ } + /* what type of output files? */ + MODE = SQUARES; + if ((N = get_argument (argc, argv, "-triangles"))) { + MODE = TRIANGLES; + remove_argument (N, &argc, argv); + } + if (argc == 2) return (TRUE); - fprintf (stderr, "USAGE: skycells (level)\n"); + fprintf (stderr, "USAGE: skycells (level) [-triangles] [-scale arcsec/pix]\n"); exit (2); } Index: /trunk/Ohana/src/addstar/src/sky_tessalation.c =================================================================== --- /trunk/Ohana/src/addstar/src/sky_tessalation.c (revision 12772) +++ /trunk/Ohana/src/addstar/src/sky_tessalation.c (revision 12773) @@ -4,9 +4,30 @@ static Coords *refcoords = NULL; -Image *sky_tessalation (int level, int *nimages) { - - int i, N, Ntriangles; +Image *sky_tessalation (int level, int *nimages, int mode) { + + Image *images; + + if (mode == SQUARES) { + images = sky_tessalation_squares (level, nimages); + return images; + } + + if (mode == TRIANGLES) { + images = sky_tessalation_triangles (level, nimages); + return images; + } + + return NULL; +} + +Image *sky_tessalation_triangles (int level, int *nimages) { + + int i, Ndigit, Ntriangles; SkyTriangle *tri, *new; - Image *image, *outimage; + Image *image; + char format[16]; + + Ndigit = (int)(log10(20*pow(4.0, level))) + 1 ; + snprintf (format, 16, "skytri.%%0%dd", Ndigit); sky_tessalation_init (); @@ -22,28 +43,46 @@ ALLOCATE (image, Image, Ntriangles); for (i = 0; i < Ntriangles; i++) { - // XXX add a user-option to choose betwen these - // sky_triangle_to_image (&image[i], &tri[i]); + sky_triangle_to_image (&image[i], &tri[i]); + snprintf (image[i].name, 32, format, i); + } + + *nimages = Ntriangles; + return (image); +} + +Image *sky_tessalation_squares (int level, int *nimages) { + + int i, N, Ndigit, Ntriangles; + SkyTriangle *tri, *new; + Image *image, *outimage; + char format[16]; + + // XXX I need to adjust for the number of cell subdivisions + Ndigit = (int)(log10(20*pow(4.0, level))) + 1 ; + snprintf (format, 16, "skycell.%%0%dd", Ndigit); + + sky_tessalation_init (); + + tri = sky_base_triangles (&Ntriangles); + + for (i = 0; i < level; i++) { + new = sky_divide_triangles (tri, &Ntriangles); + free (tri); + tri = new; + } + + ALLOCATE (image, Image, Ntriangles); + for (i = 0; i < Ntriangles; i++) { sky_triangle_to_rectangle (&image[i], &tri[i]); } - FILE *f1 = fopen ("par.1.dat", "w"); - FILE *f2 = fopen ("par.2.dat", "w"); - ALLOCATE (outimage, Image, Ntriangles); for (i = N = 0; i < Ntriangles; i++) { - if (!strcmp(image[i].coords.ctype, "DROP")) { - fprintf (f2, "%f %f\n", image[i].coords.crval1, image[i].coords.crval2); - continue; - } else { - fprintf (f1, "%f %f\n", image[i].coords.crval1, image[i].coords.crval2); - } + if (!strcmp(image[i].coords.ctype, "DROP")) continue; memcpy (&outimage[N], &image[i], sizeof(Image)); - sprintf (outimage[N].name, "test.%03d", i); + snprintf (outimage[N].name, 32, format, i); N++; } free (image); - - fclose (f1); - fclose (f2); *nimages = N; @@ -124,4 +163,5 @@ double xo, yo, xc, yc, angle, scale; double Xmin, Xmax, Ymin, Ymax; + double dB, dP, s1, s2, r1, r2, dr, dx, dy; // calculate the triangle coordinates in r,d @@ -149,4 +189,7 @@ assert (peak != -1); + // angle is from the center to the peak corner + angle = atan2(parity*xv[peak], parity*yv[peak]); // note that this is x/y not y/x (and in radians) + // find the base and height b1 = (peak + 1) % 3; @@ -156,4 +199,26 @@ xo = 0.5*(xv[b2] + xv[b1]); yo = 0.5*(yv[b2] + yv[b1]); + + // find the max perpendicular distance from the peak + + // dB[b1] == dB[b2] (since xo,yo is the midpoint of [b1] to [b2] + dB = hypot(xo -xv[b1], yo -yv[b1]); + dP = hypot(xv[peak]-xo, yv[peak] -yo); + + // XXX we could just choose the point based on s1 vs s2... + s1 = hypot(xv[peak]-xv[b1], yv[peak]-yv[b1]); + s2 = hypot(xv[peak]-xv[b2], yv[peak]-yv[b2]); + + r1 = (SQ(s1) - SQ(dB) - SQ(dP)) / (2*dP); + r2 = (SQ(s2) - SQ(dB) - SQ(dP)) / (2*dP); + + // dr >= 0 + dr = MAX (r1, r2); + + dx = -parity*dr*sin(angle); + dy = -parity*dr*cos(angle); + + xo += dx; + yo += dy; // xc, yc is the true image center @@ -163,5 +228,5 @@ NX = hypot((xv[b2]-xv[b1]),(yv[b2]-yv[b1])); NY = hypot((xv[peak]-xo),(yv[peak]-yo)); - angle = atan2(parity*xv[peak], parity*yv[peak]); // note that this is x/y not y/x (and in radians) + memset (image, 0, sizeof(Image)); Index: /trunk/Ohana/src/addstar/src/skycells.c =================================================================== --- /trunk/Ohana/src/addstar/src/skycells.c (revision 12772) +++ /trunk/Ohana/src/addstar/src/skycells.c (revision 12773) @@ -17,5 +17,5 @@ level = atoi (argv[1]); - images = sky_tessalation (level, &Nimages); + images = sky_tessalation (level, &Nimages, MODE); /*** update the image table ***/