IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Mar 8, 2007, 12:21:18 PM (19 years ago)
Author:
eugene
Message:

complete merge from branch dvo-mods-2007-02

File:
1 edited

Legend:

Unmodified
Added
Removed

  • trunk/Ohana/src/addstar/src/find_matches.c

    r10939 r12332  
    44
    55  int i, j, n, N, J, status, Nstars;
    6   double X, Y, RADIUS, RADIUS2;
     6  double RADIUS, RADIUS2;
    77  float *X1, *Y1, *X2, *Y2;
    88  float dX, dY, dR;
    … …  
    1010  int Nave, NAVE, Nmeas, NMEAS, Nmatch;
    1111  int Nsecfilt, Nsec;
    12   float Mcat, *Mval;
    13   PhotCode *code;
    1412  Coords tcoords;
    1513
    16   /* photcode data - must by of type DEP; options.photcode is equiv PRI/SEC photcode */
    17   /* this function requires incoming stars to have the same photcode.equiv value */
     14  /* photcode data - must by of type DEP; options.photcode is equiv photcode for all input
     15     images this function requires incoming stars to have the same photcode.equiv value.  if
     16     this value is not a valid photcode (ie, 0), then no modification is made to the average
     17     magnitudes (Nsec will be -1) */
     18
    1819  Nsecfilt = GetPhotcodeNsecfilt ();
    1920  Nsec     = GetPhotcodeNsec (options.photcode);
    … …  
    142143      /** add measurements for this star **/
    143144      /** dR,dD now represent arcsec **/
    144       catalog[0].measure[Nmeas].dR_PS       = 3600.0*(catalog[0].average[n].R - stars[N].R);
    145       if (catalog[0].measure[Nmeas].dR_PS > +180.0*3600.0) {
     145      catalog[0].measure[Nmeas].dR       = 3600.0*(catalog[0].average[n].R - stars[N].R);
     146      if (catalog[0].measure[Nmeas].dR > +180.0*3600.0) {
    146147          // average on high end of boundary, move star up
    147148          stars[N].R += 360.0;
    148           catalog[0].measure[Nmeas].dR_PS = 3600.0*(catalog[0].average[n].R - stars[N].R);
    149       }
    150       if (catalog[0].measure[Nmeas].dR_PS < -180.0*3600.0) {
     149          catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N].R);
     150      }
     151      if (catalog[0].measure[Nmeas].dR < -180.0*3600.0) {
    151152          // average on low end of boundary, move star down
    152153          stars[N].R -= 360.0;
    153           catalog[0].measure[Nmeas].dR_PS = 3600.0*(catalog[0].average[n].R - stars[N].R);
    154       }
    155       if (fabs(catalog[0].measure[Nmeas].dR_PS) > 10*RADIUS) {
     154          catalog[0].measure[Nmeas].dR = 3600.0*(catalog[0].average[n].R - stars[N].R);
     155      }
     156      if (fabs(catalog[0].measure[Nmeas].dR) > 10*RADIUS) {
    156157          fprintf (stderr, "error: %10.6f,%10.6f vs %10.6f,%10.6f (%f,%f vs %f,%f)\n",
    157158                   catalog[0].average[n].R, catalog[0].average[n].D,
    … …  
    160161                   Y1[i], Y2[J]);
    161162      }
    162       catalog[0].measure[Nmeas].dD_PS       = 3600.0*(catalog[0].average[n].D - stars[N].D);
    163       catalog[0].measure[Nmeas].M_PS        = stars[N].M;
    164       catalog[0].measure[Nmeas].dM_PS       = stars[N].dM;  /* error in input files stored in thousandths of mag */
    165       catalog[0].measure[Nmeas].Mcal_PS     = stars[N].Mcal;
    166       catalog[0].measure[Nmeas].t           = stars[N].t;
    167       catalog[0].measure[Nmeas].averef      = n;              /* this must be an absolute sequence number, if partial average is loaded */
    168       catalog[0].measure[Nmeas].source      = stars[N].code;  /* photcode */
    169       catalog[0].measure[Nmeas].dophot      = stars[N].dophot; 
    170       catalog[0].measure[Nmeas].flags       = 0;
    171       catalog[0].measure[Nmeas].dt_PS       = stars[N].dt;
    172       catalog[0].measure[Nmeas].airmass_PS  = stars[N].airmass;
    173 
    174       catalog[0].measure[Nmeas].Mgal_PS     = stars[N].Mgal;
    175       catalog[0].measure[Nmeas].FWx         = MIN (100*stars[N].fx, NO_MAG);
    176       catalog[0].measure[Nmeas].FWy         = MIN (100*stars[N].fy, NO_MAG);
    177       catalog[0].measure[Nmeas].theta       = MIN ((255/360)*stars[N].df, NO_ERR);
     163      catalog[0].measure[Nmeas].dD       = 3600.0*(catalog[0].average[n].D - stars[N].D);
     164
     165      /* XXX need to add dX, dY : need to load into stars[N].dX,dY */
     166      /* XXX need to add stargal, Sky, dSky, qPSF, detID, imageID */
     167      catalog[0].measure[Nmeas].Xccd     = stars[N].X;
     168      catalog[0].measure[Nmeas].Yccd     = stars[N].Y;
     169
     170      catalog[0].measure[Nmeas].M        = stars[N].M;
     171      catalog[0].measure[Nmeas].dM       = stars[N].dM;  /* error in input files stored in thousandths of mag */
     172      catalog[0].measure[Nmeas].Mcal     = stars[N].Mcal;
     173      catalog[0].measure[Nmeas].t        = stars[N].t;
     174      catalog[0].measure[Nmeas].averef   = n;              /* this must be an absolute sequence number, if partial average is loaded */
     175      catalog[0].measure[Nmeas].photcode = stars[N].code;  /* photcode */
     176      catalog[0].measure[Nmeas].dophot   = stars[N].dophot; 
     177      catalog[0].measure[Nmeas].flags    = 0;
     178      catalog[0].measure[Nmeas].dt       = stars[N].dt;
     179      catalog[0].measure[Nmeas].airmass  = stars[N].airmass;
     180
     181      catalog[0].measure[Nmeas].Mgal     = stars[N].Mgal;
     182      catalog[0].measure[Nmeas].FWx      = MIN (100*stars[N].fx, NO_MAG);
     183      catalog[0].measure[Nmeas].FWy      = MIN (100*stars[N].fy, NO_MAG);
     184      catalog[0].measure[Nmeas].theta    = MIN ((0xffff/360.0)*stars[N].df, NO_MAG);
    178185       
    179       /* it is not valid to pass PRI/SEC/REF photcodes to this routine */
    180       /* check for entries in the secfilt lists */
    181       Mcat = PhotCat (&catalog[0].measure[Nmeas]);
    182       Mval = (Nsec == -1) ? &catalog[0].average[n].M : &catalog[0].secfilt[n*Nsecfilt+Nsec].M_PS;
    183       if (*Mval == NO_MAG) *Mval = Mcat;
     186      /* set the average magnitude if not already set and if photcode.equiv is not 0 */
    184187      /* in UPDATE mode, this value is not saved; use relphot to recalculate */
     188      if (Nsec > -1) {
     189          if (catalog[0].secfilt[n*Nsecfilt+Nsec].M != NO_MAG) {
     190              catalog[0].secfilt[n*Nsecfilt+Nsec].M = PhotCat (&catalog[0].measure[Nmeas]);
     191          }
     192      }
    185193
    186194      /* adds the measurement to the calibration if appropriate color terms are found */
    … …  
    242250    catalog[0].average[Nave].R         = stars[N].R;
    243251    catalog[0].average[Nave].D         = stars[N].D;
    244     catalog[0].average[Nave].M         = NO_MAG;
    245     catalog[0].average[Nave].dM        = NO_MAG;
     252    catalog[0].average[Nave].Xp        = 0;
    246253    catalog[0].average[Nave].Nm        = 1;
    247254    catalog[0].average[Nave].Nn        = 0;
    248     catalog[0].average[Nave].Xp        = NO_MAG;
    249     catalog[0].average[Nave].Xm        = NO_MAG;
    250     catalog[0].average[Nave].Xg        = NO_MAG;
    251255    catalog[0].average[Nave].offset    = Nmeas;
    252256    catalog[0].average[Nave].missing   = -1;
    … …  
    263267
    264268    for (j = 0; j < Nsecfilt; j++) {
    265       catalog[0].secfilt[Nave*Nsecfilt+j].M_PS  = NO_MAG;
    266       catalog[0].secfilt[Nave*Nsecfilt+j].dM_PS = NO_MAG;
    267       catalog[0].secfilt[Nave*Nsecfilt+j].Xm    = NO_MAG;
    268     }
    269 
    270     catalog[0].measure[Nmeas].dR_PS       = 0.0;
    271     catalog[0].measure[Nmeas].dD_PS       = 0.0;
    272     catalog[0].measure[Nmeas].M_PS        = stars[N].M;
    273     catalog[0].measure[Nmeas].dM_PS       = stars[N].dM;
    274     catalog[0].measure[Nmeas].Mcal_PS     = stars[N].Mcal;
    275     catalog[0].measure[Nmeas].t           = stars[N].t;
    276     catalog[0].measure[Nmeas].averef      = Nave;           /* XXX EAM : must be absolute Nave if partial read */
    277     catalog[0].measure[Nmeas].source      = stars[N].code;  /* photcode */
    278     catalog[0].measure[Nmeas].dophot      = stars[N].dophot; 
    279     catalog[0].measure[Nmeas].flags       = 0;
    280     catalog[0].measure[Nmeas].dt_PS       = stars[N].dt;
    281     catalog[0].measure[Nmeas].airmass_PS  = stars[N].airmass;
    282 
    283     catalog[0].measure[Nmeas].Mgal_PS     = stars[N].Mgal;
    284     catalog[0].measure[Nmeas].FWx         = MIN (100*stars[N].fx, NO_MAG);
    285     catalog[0].measure[Nmeas].FWy         = MIN (100*stars[N].fy, NO_MAG);
    286     catalog[0].measure[Nmeas].theta       = MIN ((255/360)*stars[N].df, NO_ERR);
    287 
    288     Mcat = PhotCat (&catalog[0].measure[Nmeas]);
    289     Mval = (Nsec == -1) ? &catalog[0].average[Nave].M : &catalog[0].secfilt[Nave*Nsecfilt+Nsec].M_PS;
    290     if (*Mval == NO_MAG) *Mval = Mcat;
     269      catalog[0].secfilt[Nave*Nsecfilt+j].M  = NO_MAG;
     270      catalog[0].secfilt[Nave*Nsecfilt+j].dM = NO_MAG;
     271      catalog[0].secfilt[Nave*Nsecfilt+j].Xm = NO_MAG;
     272    }
     273
     274    /* XXX need to add dX, dY : need to load into stars[N].dX,dY */
     275    /* XXX need to add stargal, Sky, dSky, qPSF, detID, imageID */
     276    catalog[0].measure[Nmeas].Xccd     = stars[N].X;
     277    catalog[0].measure[Nmeas].Yccd     = stars[N].Y;
     278
     279    catalog[0].measure[Nmeas].dR        = 0.0;
     280    catalog[0].measure[Nmeas].dD        = 0.0;
     281    catalog[0].measure[Nmeas].M         = stars[N].M;
     282    catalog[0].measure[Nmeas].dM        = stars[N].dM;
     283    catalog[0].measure[Nmeas].Mcal      = stars[N].Mcal;
     284    catalog[0].measure[Nmeas].t         = stars[N].t;
     285    catalog[0].measure[Nmeas].averef    = Nave;           /* XXX EAM : must be absolute Nave if partial read */
     286    catalog[0].measure[Nmeas].photcode  = stars[N].code;  /* photcode */
     287    catalog[0].measure[Nmeas].dophot    = stars[N].dophot; 
     288    catalog[0].measure[Nmeas].flags     = 0;
     289    catalog[0].measure[Nmeas].dt        = stars[N].dt;
     290    catalog[0].measure[Nmeas].airmass   = stars[N].airmass;
     291
     292    catalog[0].measure[Nmeas].Mgal      = stars[N].Mgal;
     293    catalog[0].measure[Nmeas].FWx       = MIN (100*stars[N].fx, NO_MAG);
     294    catalog[0].measure[Nmeas].FWy       = MIN (100*stars[N].fy, NO_MAG);
     295    catalog[0].measure[Nmeas].theta     = MIN ((0xffff/360.0)*stars[N].df, NO_MAG);
     296    /* XXX replace df here with theta, right? */
     297
     298    /* set the average magnitude if not already set and the photcode.equiv is not 0 */
     299    /* in UPDATE mode, this value is not saved; use relphot to recalculate */
     300    if (Nsec > -1) {
     301        catalog[0].secfilt[Nave*Nsecfilt+Nsec].M = PhotCat (&catalog[0].measure[Nmeas]);
     302    }
    291303
    292304    /* next[Nmeas] should always be -1 in this context (it is always the only
Note: See TracChangeset for help on using the changeset viewer.