IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 12751


Ignore:
Timestamp:
Apr 5, 2007, 9:43:25 AM (19 years ago)
Author:
eugene
Message:

fixed bugs

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/Ohana/src/addstar/src/sky_tessalation.c

    r12749 r12751  
    3232int sky_triangle_to_image (Image *image, SkyTriangle *triangle) {
    3333
    34   int i, peak, b1, b2, NX, NY;
     34  int i, n, parity, peak, b1, b2, NX, NY;
    3535  double xv[3], yv[3];        // coordinates of the vertex in the reference projection
    36   double xo, yo, angle, scale;
     36  double xo, yo, xc, yc, angle, scale;
    3737
    3838  // calculate the triangle coordinates in r,d
    3939  sky_triangle_coords (triangle);
    4040
    41   // we will project to the triangle center position
    42   refcoords[0].crval1 = triangle[0].r;
    43   refcoords[0].crval2 = triangle[0].d;
    44 
    45   // find the size, rotation, and parity of the image
    46   // project the vertices and find the image parity
    47   int parity = 1;
    48   for (i = 0; i < 3; i++) {
    49     RD_to_XY (&xv[i], &yv[i], triangle[0].rv[i], triangle[0].dv[i], refcoords);
    50     parity *= SIGN(yv[i]);
    51   }
    52 
    53   // choose the peak vertex
    54   peak = -1;
    55   for (i = 0; (peak == -1) && (i < 3); i++) {
    56     if (parity == SIGN(yv[i])) {
    57       peak = i;
    58     }
    59   }
    60   assert (peak != -1);
    61   angle = atan2(xv[peak], yv[peak]); // note that this is x/y not y/x (and in radians)
    62 
    63   // find the base and height
    64   b1 = (peak + 1) % 3;
    65   b2 = (peak + 2) % 3;
     41  for (n = 0; n < 3; n++) {
     42    // we will project to the triangle center position
     43    refcoords[0].crval1 = triangle[0].r;
     44    refcoords[0].crval2 = triangle[0].d;
     45
     46    // find the size, rotation, and parity of the image
     47    // project the vertices and find the image parity
     48    parity = 1;
     49    for (i = 0; i < 3; i++) {
     50      RD_to_XY (&xv[i], &yv[i], triangle[0].rv[i], triangle[0].dv[i], refcoords);
     51      parity *= SIGN(yv[i]);
     52    }
     53
     54    // choose the peak vertex
     55    peak = -1;
     56    for (i = 0; (peak == -1) && (i < 3); i++) {
     57      if (parity == SIGN(yv[i])) {
     58        peak = i;
     59      }
     60    }
     61    assert (peak != -1);
     62
     63    // find the base and height
     64    b1 = (peak + 1) % 3;
     65    b2 = (peak + 2) % 3;
     66
     67    // xo, yo is the center of the baseline
     68    xo = 0.5*(xv[b2] + xv[b1]);
     69    yo = 0.5*(yv[b2] + yv[b1]);
     70
     71    if (n != 2) {
     72      // xc, yc is the true image center
     73      xc = 0.5*(xv[peak] + xo);
     74      yc = 0.5*(yv[peak] + yo);
     75      XY_to_RD (&triangle[0].r, &triangle[0].d, xc, yc, refcoords);
     76    }
     77  } 
    6678
    6779  NX = hypot((xv[b2]-xv[b1]),(yv[b2]-yv[b1]));
    68   xo = 0.5*(xv[b2] + xv[b1]);
    69   yo = 0.5*(yv[b2] + yv[b1]);
    7080  NY = hypot((xv[peak]-xo),(yv[peak]-yo));
     81  angle = atan2(parity*xv[peak], parity*yv[peak]); // note that this is x/y not y/x (and in radians)
    7182
    7283  memset (image, 0, sizeof(Image));
    7384  image[0].coords = *refcoords;
    7485  image[0].coords.pc1_1 = +cos(angle);
    75   image[0].coords.pc1_2 = -sin(angle);
    76   image[0].coords.pc2_1 = +sin(angle);
     86  image[0].coords.pc1_2 = +sin(angle);
     87  image[0].coords.pc2_1 = -sin(angle);
    7788  image[0].coords.pc2_2 = +cos(angle);
    7889
     
    8394  }
    8495  sprintf (image[0].name, "test.%03d", i);
    85   fprintf (stderr, "NX, NY: %d, %d\n", NX, NY);
    8696  if ((NX > 60000) || (NY > 60000)) {
    8797    scale = NX / 60000.0;
     
    95105  image[0].NY = NY;
    96106  image[0].code = 1;
     107
     108  // fprintf (stderr, "%d  %f  %d %d %d  %d %d\n", parity, angle*DEG_RAD, peak, b1, b2, NX, NY);
     109
    97110  return (TRUE);
    98111}
     
    122135    triangle[0].rv[i] = DEG_RAD * atan2(triangle[0].vertex[i].y, triangle[0].vertex[i].x);
    123136  }
     137
     138  return TRUE;
    124139}
    125140
     
    128143SkyTriangle *sky_divide_triangles (SkyTriangle *in, int *ntriangles) {
    129144
    130   int i, j, Ntriangles;
     145  int i, j, Ntriangles, Nt;
    131146  SkyTriangle *out;
    132147
     
    134149  ALLOCATE (out, SkyTriangle, Ntriangles);
    135150
     151  Nt = 0;
    136152  for (i = 0; i < *ntriangles; i++) {
    137153    for (j = 0; j < 3; j++) {
     
    145161      out[4*i + 3].vertex[j] = out[4*i + j].vertex[1];
    146162    }
    147   }
    148   *ntriangles = Ntriangles;
     163# if (0)
     164    fprintf (stdout, "%f %f %f  %f %f %f  %f %f %f\n",
     165             in[i].vertex[0].x, in[i].vertex[0].y, in[i].vertex[0].z,
     166             in[i].vertex[1].x, in[i].vertex[1].y, in[i].vertex[1].z,
     167             in[i].vertex[2].x, in[i].vertex[2].y, in[i].vertex[2].z);
     168    for (j = 0; j < 4; j++) {
     169      fprintf (stdout, "%f %f %f  %f %f %f  %f %f %f\n",
     170               out[4*i+j].vertex[0].x, out[4*i+j].vertex[0].y, out[4*i+j].vertex[0].z,
     171               out[4*i+j].vertex[1].x, out[4*i+j].vertex[1].y, out[4*i+j].vertex[1].z,
     172               out[4*i+j].vertex[2].x, out[4*i+j].vertex[2].y, out[4*i+j].vertex[2].z);
     173    }
     174# endif 
     175    Nt += 4;
     176  }
     177  *ntriangles = Nt;
    149178  return (out);
    150179}
     
    239268    tri[i].vertex[1].x = +ctht*cos(psi);
    240269    tri[i].vertex[1].y = +ctht*sin(psi);
    241     tri[i].vertex[1].z = +stht;
     270    tri[i].vertex[1].z = -stht;
    242271   
    243272    psi = (i + 1.5)*D_PSI;
    244273    tri[i].vertex[2].x = +ctht*cos(psi);
    245274    tri[i].vertex[2].y = +ctht*sin(psi);
    246     tri[i].vertex[2].z = +stht;
     275    tri[i].vertex[2].z = -stht;
    247276
    248277    tri[i].vertex[0].x = +0;
Note: See TracChangeset for help on using the changeset viewer.