IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Mar 6, 2011, 2:34:55 PM (15 years ago)
Author:
eugene
Message:

updates from trunk

Location:
branches/eam_branches/ipp-20110213/Ohana/src
Files:
4 edited

Legend:

Unmodified
Added
Removed
  • branches/eam_branches/ipp-20110213/Ohana/src/addstar/src/ReadStarsFITS.c

    r30613 r30812  
    4545  }
    4646  if (!strcmp (type, "PS1_V2")) {
     47    if (table.header[0].Naxis[0] == 20) {
     48      // skip the invalid DETEFF tables which were mistakenly labeled as PS1_V2
     49      return (NULL);
     50    }
    4751    stars = Convert_PS1_V2 (&table, &Nstars);
    4852  }
     
    5155  }
    5256  if (stars == NULL) {
    53     fprintf (stderr, "ERROR: invalid table type %s\n", type);
     57    fprintf (stderr, "invalid table type %s\n", type);
    5458    return (NULL);
    5559  }
     
    212216  Stars *stars;
    213217  CMF_PS1_V1 *ps1data;
     218
     219  // CMF_PS1_V1 was modified 2009.05.26 (r24251) to use doubles for ra & dec.  this was a
     220  // mistake in two ways: a new format should have been defined (eg, CMF_PS1_V2), and the
     221  // layout used did not have clean byte-boundaries for the corresponding structure.  The
     222  // former means we have two varieties of CMF_PS1_V1 out there; the latter means that the
     223  // autocode tools do not work to read in the new version, even if we recognize it.  Here we
     224  // test for the existence of the broken version (table[0].headers[0].Naxis[0] == 136), and
     225  // call a special conversion function if it is found.
     226
     227  if (table[0].header[0].Naxis[0] == 136) {
     228      stars = Convert_PS1_V1_Alt (table, nstars);
     229      return (stars);
     230  }
     231
     232  ps1data = gfits_table_get_CMF_PS1_V1 (table, &Nstars, NULL);
     233  if (!ps1data) {
     234    fprintf (stderr, "skipping inconsistent entry\n");
     235    return (NULL);
     236  }
     237  ZeroPt = GetZeroPoint();
     238
     239  ALLOCATE (stars, Stars, Nstars);
     240  for (i = 0; i < Nstars; i++) {
     241    InitStar (&stars[i]);
     242    stars[i].measure.Xccd       = ps1data[i].X;
     243    stars[i].measure.Yccd       = ps1data[i].Y;
     244    stars[i].measure.dXccd      = ToShortPixels(ps1data[i].dX);
     245    stars[i].measure.dYccd      = ToShortPixels(ps1data[i].dY);
     246
     247    stars[i].measure.posangle   = ToShortDegrees(ps1data[i].posangle);
     248    stars[i].measure.pltscale   = ps1data[i].pltscale;
     249
     250    if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
     251        stars[i].measure.M      = NAN;
     252    } else {
     253        stars[i].measure.M      = ps1data[i].M + ZeroPt;
     254    }
     255    stars[i].measure.dM         = ps1data[i].dM;
     256    stars[i].measure.dMcal      = ps1data[i].dMcal;
     257    stars[i].measure.Map        = ps1data[i].Map + ZeroPt;
     258                       
     259    stars[i].measure.Sky        = ps1data[i].sky;
     260    stars[i].measure.dSky       = ps1data[i].dSky;
     261                       
     262    stars[i].measure.psfChisq   = ps1data[i].psfChisq;
     263    stars[i].measure.psfQual    = ps1data[i].psfQual;
     264    stars[i].measure.psfNdof    = ps1data[i].psfNdof;
     265    stars[i].measure.psfNpix    = ps1data[i].psfNpix;
     266    stars[i].measure.crNsigma   = ps1data[i].crNsigma;
     267    stars[i].measure.extNsigma  = ps1data[i].extNsigma;
     268
     269    stars[i].measure.FWx        = ToShortPixels(ps1data[i].fx);
     270    stars[i].measure.FWy        = ToShortPixels(ps1data[i].fy);
     271    stars[i].measure.theta      = ToShortDegrees(ps1data[i].df);
     272
     273    stars[i].measure.Mxx        = ToShortPixels(ps1data[i].Mxx);
     274    stars[i].measure.Mxy        = ToShortPixels(ps1data[i].Mxy);
     275    stars[i].measure.Myy        = ToShortPixels(ps1data[i].Myy);
     276                       
     277    stars[i].measure.photFlags  = ps1data[i].flags;
     278
     279    // this is may optionally be replaced by the internal sequence (see FilterStars.c)
     280    stars[i].measure.detID      = ps1data[i].detID;
     281
     282    // the Average fields and the following Measure fields are set in FilterStars after
     283    // the image metadata is in hand:  dR, dD, Mcal, dt, airmass, az, t, imageID, extID,
     284    // averef is set in find_matches, dbFlags is zero on ingest.
     285
     286    // the following fields are currently not being set anywhere: t_msec
     287  }   
     288  *nstars = Nstars;
     289  return (stars);
     290}
     291
     292Stars *Convert_PS1_V1_Alt (FTable *table, unsigned int *nstars) {
     293
     294  off_t Nstars;
     295  unsigned int i;
     296  double ZeroPt;
     297  Stars *stars;
     298  CMF_PS1_V2 *ps1data;
    214299
    215300  // CMF_PS1_V1 was modified 2009.05.26 (r24251) to use doubles for ra & dec.  this was a
     
    221306  // call a special conversion function if it is found.
    222307
    223   if (table[0].header[0].Naxis[0] == 136) {
    224       stars = Convert_PS1_V1_Alt (table, nstars);
    225       return (stars);
    226   }
    227 
    228   ps1data = gfits_table_get_CMF_PS1_V1 (table, &Nstars, NULL);
     308  ps1data = gfits_table_get_CMF_PS1_V1_Alt (table, &Nstars, NULL);
    229309  if (!ps1data) {
    230310    fprintf (stderr, "skipping inconsistent entry\n");
     
    286366}
    287367
    288 Stars *Convert_PS1_V1_Alt (FTable *table, unsigned int *nstars) {
     368Stars *Convert_PS1_V2 (FTable *table, unsigned int *nstars) {
    289369
    290370  off_t Nstars;
     
    294374  CMF_PS1_V2 *ps1data;
    295375
    296   // CMF_PS1_V1 was modified 2009.05.26 (r24251) to use doubles for ra & dec.  this was a
    297   // mistake in two ways: a few format should have been defined (eg, CMF_PS1_V2), and the
    298   // layout used did not have clean byte-boundaries for the corresponding structure.  The
    299   // former means we have two varieties of CMF_PS1_V1 out there; the latter means that the
    300   // autocode tools do not work to read in the new version, even if we recognize it.  Here we
    301   // test for the existence of the broken version (table[0].headers[0].Naxis[0] == 136), and
    302   // call a special conversion function if it is found.
    303 
    304   ps1data = gfits_table_get_CMF_PS1_V1_Alt (table, &Nstars, NULL);
     376  ps1data = gfits_table_get_CMF_PS1_V2 (table, &Nstars, NULL);
    305377  if (!ps1data) {
    306378    fprintf (stderr, "skipping inconsistent entry\n");
     
    362434}
    363435
    364 Stars *Convert_PS1_V2 (FTable *table, unsigned int *nstars) {
     436Stars *Convert_PS1_V3 (FTable *table, unsigned int *nstars) {
    365437
    366438  off_t Nstars;
     
    368440  double ZeroPt;
    369441  Stars *stars;
    370   CMF_PS1_V2 *ps1data;
    371 
    372   ps1data = gfits_table_get_CMF_PS1_V2 (table, &Nstars, NULL);
     442  CMF_PS1_V3 *ps1data;
     443
     444  ps1data = gfits_table_get_CMF_PS1_V3 (table, &Nstars, NULL);
    373445  if (!ps1data) {
    374446    fprintf (stderr, "skipping inconsistent entry\n");
     
    429501  return (stars);
    430502}
    431 
    432 Stars *Convert_PS1_V3 (FTable *table, unsigned int *nstars) {
    433 
    434   off_t Nstars;
    435   unsigned int i;
    436   double ZeroPt;
    437   Stars *stars;
    438   CMF_PS1_V3 *ps1data;
    439 
    440   ps1data = gfits_table_get_CMF_PS1_V3 (table, &Nstars, NULL);
    441   if (!ps1data) {
    442     fprintf (stderr, "skipping inconsistent entry\n");
    443     return (NULL);
    444   }
    445   ZeroPt = GetZeroPoint();
    446 
    447   ALLOCATE (stars, Stars, Nstars);
    448   for (i = 0; i < Nstars; i++) {
    449     InitStar (&stars[i]);
    450     stars[i].measure.Xccd       = ps1data[i].X;
    451     stars[i].measure.Yccd       = ps1data[i].Y;
    452     stars[i].measure.dXccd      = ToShortPixels(ps1data[i].dX);
    453     stars[i].measure.dYccd      = ToShortPixels(ps1data[i].dY);
    454 
    455     stars[i].measure.posangle   = ToShortDegrees(ps1data[i].posangle);
    456     stars[i].measure.pltscale   = ps1data[i].pltscale;
    457 
    458     if ((ps1data[i].M >= 0.0) || isnan(ps1data[i].M)) {
    459         stars[i].measure.M      = NAN;
    460     } else {
    461         stars[i].measure.M      = ps1data[i].M + ZeroPt;
    462     }
    463     stars[i].measure.dM         = ps1data[i].dM;
    464     stars[i].measure.dMcal      = ps1data[i].dMcal;
    465     stars[i].measure.Map        = ps1data[i].Map + ZeroPt;
    466                        
    467     stars[i].measure.Sky        = ps1data[i].sky;
    468     stars[i].measure.dSky       = ps1data[i].dSky;
    469                        
    470     stars[i].measure.psfChisq   = ps1data[i].psfChisq;
    471     stars[i].measure.psfQual    = ps1data[i].psfQual;
    472     stars[i].measure.psfNdof    = ps1data[i].psfNdof;
    473     stars[i].measure.psfNpix    = ps1data[i].psfNpix;
    474     stars[i].measure.crNsigma   = ps1data[i].crNsigma;
    475     stars[i].measure.extNsigma  = ps1data[i].extNsigma;
    476 
    477     stars[i].measure.FWx        = ToShortPixels(ps1data[i].fx);
    478     stars[i].measure.FWy        = ToShortPixels(ps1data[i].fy);
    479     stars[i].measure.theta      = ToShortDegrees(ps1data[i].df);
    480 
    481     stars[i].measure.Mxx        = ToShortPixels(ps1data[i].Mxx);
    482     stars[i].measure.Mxy        = ToShortPixels(ps1data[i].Mxy);
    483     stars[i].measure.Myy        = ToShortPixels(ps1data[i].Myy);
    484                        
    485     stars[i].measure.photFlags  = ps1data[i].flags;
    486 
    487     // this is may optionally be replaced by the internal sequence (see FilterStars.c)
    488     stars[i].measure.detID      = ps1data[i].detID;
    489 
    490     // the Average fields and the following Measure fields are set in FilterStars after
    491     // the image metadata is in hand:  dR, dD, Mcal, dt, airmass, az, t, imageID, extID,
    492     // averef is set in find_matches, dbFlags is zero on ingest.
    493 
    494     // the following fields are currently not being set anywhere: t_msec
    495   }   
    496   *nstars = Nstars;
    497   return (stars);
    498 }
  • branches/eam_branches/ipp-20110213/Ohana/src/dvomerge/src/merge_catalogs_old.c

    r29938 r30812  
    219219      }
    220220      if (fabs(output[0].measure[Nmeas].dR) > 10*RADIUS) {
    221         fprintf (stderr, "error: %10.6f,%10.6f vs %10.6f,%10.6f (%f,%f vs %f,%f)\n",
    222                  output[0].average[n].R, output[0].average[n].D, Rin, Din,
    223                  X1[i], X2[Jmin], Y1[i], Y2[Jmin]);
    224         // XXX abort on this? -- this is a bad failure...
     221        // ok take declination into account and check again.
     222        double cosD = cos(RAD_DEG*Din);
     223        if (fabs(output[0].measure[Nmeas].dR*cosD) > 10*RADIUS) {
     224
     225            fprintf (stderr, "error: %10.6f,%10.6f vs %10.6f,%10.6f (%f,%f vs %f,%f)\n",
     226             output[0].average[n].R, output[0].average[n].D, Rin, Din,
     227             X1[i], X2[Jmin], Y1[i], Y2[Jmin]);
     228
     229            // XXX abort on this? -- this is a bad failure...
     230        }
    225231      }
    226232      input[0].found[N] = Nmeas;
  • branches/eam_branches/ipp-20110213/Ohana/src/getstar/src/MatchCoords.c

    r30499 r30812  
    33# define FLT_MAX 1e32
    44# endif
     5int ComputeLMLimits(double xmin, double xmax, double ymin, double ymax,
     6            double *Lmin, double *Lmax, double *Mmin, double *Mmax, Coords *coords);
    57
    68/* given coordinate, find images in list that contain the point */
     
    1214  double Xi[4], Yi[4];  /* image and original corners */
    1315  double xmin, xmax, ymin, ymax;
     16  double Lmin, Lmax, Mmin, Mmax;
    1417
    1518
     
    2023     also define the vtable entries for the images we keep  */
    2124
     25  int haveLimits = 0;
    2226  for (i = 0; i < NdbImages; i++) {
     27      haveLimits = 0;
    2328
    2429    if (!WITH_PHU && !strcmp (&dbImages[i].coords.ctype[4], "-DIS")) continue;
     
    3035    /* define image corners */
    3136    SetImageCorners (Xi, Yi, &dbImages[i]);
    32     // Xi[4] = Xi[0]; Yi[4] = Yi[0];
     37    // Xi[3] = Xi[0]; Yi[3] = Yi[0];
    3338
    3439    ymin = xmin = +FLT_MAX;
     
    5257            continue;
    5358        }
    54 
    5559       
    56         if (VERBOSE) fprintf(stderr, "%ld %8.1f %8.1f %s\n", i, x, y, dbImages[i].name);
     60        if (VERBOSE) fprintf(stderr, OFF_T_FMT" %8.1f %8.1f %s\n", i, x, y, dbImages[i].name);
    5761
    5862        if ((x >= xmin && x <= xmax) && (y >= ymin && y <= ymax)) {
     63            if (VERBOSE) fprintf(stderr, "\tpotential overlap for point %d image %i\n", j, (int) i);
     64
     65            // We may be using the astrometry far outside it's range of validity.
     66            // As a sanity check make sure that the focal plane coords are in bounds.
     67            // See ticket 1381
     68            if (!haveLimits) {
     69                // We defer computing the focal plane limits until we need them
     70                status = ComputeLMLimits(xmin, xmax, ymin, ymax,
     71                    &Lmin, &Lmax, &Mmin, &Mmax, &dbImages[i].coords);
     72                if (!status) {
     73                    continue;
     74                }
     75                haveLimits = 1;
     76            }
     77
     78            double L, M;
     79            status = RD_to_LM(&L, &M, pt->ra, pt->dec, &dbImages[i].coords);
     80            if (!status) {
     81                // This won't happen because RD_to_XY succeeded above
     82                continue;
     83            }
     84
     85            if (VERBOSE) fprintf(stderr, "L: %f M: %f\n", L, M);
     86
     87            if (L < Lmin || L > Lmax || M < Mmin || M > Mmax) {
     88                if (VERBOSE) fprintf(stderr, "\trejecting match for because focal plane coordinates are out of range\n");
     89                continue;
     90            }
     91            if (VERBOSE) fprintf(stderr, "\tmatch confirmed\n");
     92
    5993            totalMatches++;
    6094
     
    88122  }
    89123}
     124
     125// compute bounding box of this image in focal plane coordinates
     126int ComputeLMLimits(double xmin, double xmax, double ymin, double ymax,
     127            double *Lmin, double *Lmax, double *Mmin, double *Mmax, Coords *coords)
     128{
     129    double L0, L1, M0, M1;
     130    int status = XY_to_LM(&L0, &M0, xmin, ymin, coords);
     131    if (!status) {
     132        fprintf(stderr, "failed to transform image corner to LM\n");
     133        return 0;
     134    }
     135    status = XY_to_LM(&L1, &M1, xmax, ymax, coords);
     136    if (!status) {
     137        fprintf(stderr, "failed to transform image corner to LM\n");
     138        return 0;
     139    }
     140
     141    if (L1 > L0) {
     142        *Lmin = L0;
     143        *Lmax = L1;
     144    } else {
     145        *Lmin = L1;
     146        *Lmax = L0;
     147    }
     148    if (M1 > M0) {
     149        *Mmin = M0;
     150        *Mmax = M1;
     151    } else {
     152        *Mmin = M1;
     153        *Mmax = M0;
     154    }
     155    if (VERBOSE) fprintf(stderr, "Lmin: %f Lmax: %f Mmin: %f Mmax: %f\n", *Lmin, *Lmax, *Mmin, *Mmax);
     156    return 1;
     157}
  • branches/eam_branches/ipp-20110213/Ohana/src/libohana/include/ohana_allocate.h

    r30602 r30812  
    4949
    5050# define REALLOCATE(PTR,TYPE,SIZE) {                                    \
    51   PTR = (TYPE *) realloc(PTR,(unsigned)(MAX(((SIZE)*((int)sizeof(TYPE))),1))); \
     51  PTR = (TYPE *) realloc(PTR,(size_t)(MAX(((SIZE)*((int)sizeof(TYPE))),1))); \
    5252  if (PTR == NULL) {                                                    \
    5353    fprintf(stderr,"failed realloc at %d in %s\n", __LINE__, __FILE__); \
     
    5757  if ((NCURR) >= (SIZE)) {                              \
    5858    SIZE += DELTA;                                                      \
    59     PTR = (TYPE *) realloc(PTR,(unsigned)(MAX(((SIZE)*((int)sizeof(TYPE))),1))); \
     59    PTR = (TYPE *) realloc(PTR,(size_t)(MAX(((SIZE)*((int)sizeof(TYPE))),1))); \
    6060    if (PTR == NULL) {                                                  \
    6161      fprintf(stderr,"failed realloc increment at %d in %s\n", __LINE__, __FILE__); \
Note: See TracChangeset for help on using the changeset viewer.