IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Apr 19, 2013, 4:31:05 PM (13 years ago)
Author:
eugene
Message:

merge changes from eam_branches/ipp-20130307 : parallelize dvosecfilt, use dvo_average,measure,secfilt_init functions more broadly; move image calibration to relphot_images; add synthetic photometry referece (w-band); put relphot MARKTIME & INITTIME in macros; in dvo_clients, mextract fields which need image data use a subset image metadata table (instead of loading full Images.dat table); create FIELD and MOSAIC fields for mextract; somewhat better rules for photcode:ave & similar selections in mextract; add hpm_* concept in relastro (high-speed proper motions); fix precision errors in dvopsps; check for dvopsps exit conditions; error-bar clipping to limits of plotting window; parallelize delstar -dup-images; do NOT delete parents (because parent IDs are broken after dvomerge); minor handshake when setting up KAPA connection; avextract and related do not need to launch remote jobs on all clients if region does not include tables from all hosts; add fitplx and xsection functions to mana; enable addstar of diff cmfs (at least PS1_DV3)

Location:
trunk/Ohana
Files:
11 edited
2 copied

Legend:

Unmodified
Added
Removed

  • trunk/Ohana

  • trunk/Ohana/src/libdvo/Makefile

    r35263 r35416  
    3636$(DESTINC)/ps1_v3_defs.h \
    3737$(DESTINC)/ps1_v4_defs.h \
    38 $(DESTINC)/ps1_ref_defs.h
     38$(DESTINC)/ps1_ref_defs.h \
     39$(DESTINC)/cmf-ps1-dv3.h
    3940
    4041INCS = $(DEFS) $(DESTINC)/dvo.h $(DESTINC)/autocode.h $(DESTINC)/dvo_util.h $(DESTINC)/dvodb.h $(DESTINC)/libdvo_astro.h $(DESTINC)/convert.h $(DESTINC)/get_graphdata.h
    … …  
    8485$(SRC)/cmf-ps1-v1-alt.$(ARCH).o \
    8586$(SRC)/cmf-ps1-sv1-alt.$(ARCH).o \
     87$(SRC)/cmf-ps1-dv3.$(ARCH).o \
    8688$(SRC)/dvo_util.$(ARCH).o \
    8789$(SRC)/dbBooleanCond.$(ARCH).o          \

  • trunk/Ohana/src/libdvo/include/dvo.h

    r35102 r35416  
    129129  ID_MEAS_BLEND_MEAS_X   = 0x00001000,  // detection is within radius of multiple objects across catalogs                   
    130130  ID_MEAS_ARTIFACT       = 0x00002000,  // detection is thought to be non-astronomical                               
    131   ID_MEAS_UNDEF_5        = 0x00004000,  // unused
     131  ID_MEAS_SYNTH_MAG      = 0x00004000,  // magnitude is synthetic
    132132  ID_MEAS_PHOTOM_UBERCAL = 0x00008000,  // externally-supplied zero point from ubercal analysis
    133133  ID_MEAS_STACK_PRIMARY  = 0x00010000,  // this stack measurement is in the primary skycell
    … …  
    288288CMF_PS1_SV1 *gfits_table_get_CMF_PS1_SV1_Alt (FTable *ftable, off_t *Ndata, char *swapped);
    289289
     290// another special case : does not match byte-boundaries
     291# include "cmf-ps1-dv3.h"
     292
    290293typedef struct {
    291294  int Ncode;                                      // number of photcodes
    … …  
    717720
    718721void dvo_average_init (Average *average);
     722void dvo_averageT_init (AverageTiny *average);
    719723void dvo_secfilt_init (SecFilt *secfilt);
    720724void dvo_measure_init (Measure *measure);
     725void dvo_measureT_init (MeasureTiny *measure);
    721726
    722727# endif // DVO_H

  • trunk/Ohana/src/libdvo/include/dvodb.h

    r35263 r35416  
    125125      MEAS_YFIX,
    126126      MEAS_POS_SYS_ERR,
     127      MEAS_XFIELD,
     128      MEAS_YFIELD,
    127129      MEAS_XMOSAIC,
    128130      MEAS_YMOSAIC,
    … …  
    301303  double crval1;
    302304  double crval2;
     305  float theta;
    303306  unsigned int imageID;
    304307  unsigned int externID;
    … …  
    382385ImageMetadata *MatchImageMetadataDVO (unsigned int imageID);
    383386Coords *MatchMosaicMetadata (unsigned int imageID);
     387Coords *MatchFieldMetadata (unsigned int imageID);
    384388off_t match_image_by_ID (ImageMetadata *image, off_t Nimage, unsigned int ID);
    385389void sort_image_metadata (ImageMetadata *image, off_t Nimage);

  • trunk/Ohana/src/libdvo/src/ImageMetadata.c

    r35263 r35416  
    6363  GET_COLUMN (crval1,   "CRVAL1",         double);
    6464  GET_COLUMN (crval2,   "CRVAL2",         double);
     65  GET_COLUMN (theta,    "THETA",          float);
    6566  GET_COLUMN (Mcal,     "MCAL",           float);
    6667  GET_COLUMN (secz,     "SECZ",           float);
    … …  
    7576    image[i].crval1   = crval1[i]  ;
    7677    image[i].crval2   = crval2[i]  ;
     78    image[i].theta    = theta[i]  ;
    7779    image[i].Mcal     = Mcal[i]    ;
    7880    image[i].secz     = secz[i]    ;
    … …  
    8789  free (crval1);
    8890  free (crval2);
     91  free (theta);
    8992  free (Mcal);
    9093  free (secz);
    … …  
    112115  FTable ftable;
    113116
     117  BuildChipMatch (image, Nimage);
     118
    114119  gfits_init_header (&header);
    115120  header.extend = TRUE;
    … …  
    124129  gfits_define_bintable_column (&theader, "D", "CRVAL1", "ra at center", "degrees", 1.0, 0.0);
    125130  gfits_define_bintable_column (&theader, "D", "CRVAL2", "dec at center", "degrees", 1.0, 0.0);
     131  gfits_define_bintable_column (&theader, "E", "THETA", "camera rot angle", "degrees", 1.0, 0.0);
    126132  gfits_define_bintable_column (&theader, "E", "MCAL", "zero point offset", "magnitudes", 1.0, 0.0);
    127133  gfits_define_bintable_column (&theader, "E", "SECZ", "airmass", "none", 1.0, 0.0);
    … …  
    134140  unsigned int *imageID, *externID, *expname;
    135141  double *crval1, *crval2;
    136   float *Mcal, *Xcenter, *Ycenter, *secz;
     142  float *Mcal, *Xcenter, *Ycenter, *secz, *theta;
    137143
    138144  // create intermediate storage arrays
    … …  
    142148  ALLOCATE (crval1,   double,         Nimage);
    143149  ALLOCATE (crval2,   double,         Nimage);
     150  ALLOCATE (theta,    float,          Nimage);
    144151  ALLOCATE (Mcal,     float,          Nimage);
    145152  ALLOCATE (secz,     float,          Nimage);
    … …  
    149156  // assign the storage arrays
    150157  for (i = 0; i < Nimage; i++) {
     158    int Nmosaic = FindMosaicForImage (image, Nimage, i);
     159    if (!Nmosaic) continue;
     160    Nmosaic --;
    151161    imageID[i]  = image[i].imageID;
    152162    externID[i] = image[i].externID;
    153     crval1[i]   = image[i].coords.crval1;
    154     crval2[i]   = image[i].coords.crval2;
     163    crval1[i]   = image[Nmosaic].coords.crval1;
     164    crval2[i]   = image[Nmosaic].coords.crval2;
     165
     166    theta[i]    = DEG_RAD*atan2(image[Nmosaic].coords.pc1_2, image[Nmosaic].coords.pc1_1);
     167
    155168    Mcal[i]     = image[i].Mcal;
    156169    secz[i]     = image[i].secz;
    … …  
    172185  gfits_set_bintable_column (&theader, &ftable, "CRVAL1",         crval1,  Nimage);
    173186  gfits_set_bintable_column (&theader, &ftable, "CRVAL2",         crval2,  Nimage);
     187  gfits_set_bintable_column (&theader, &ftable, "THETA",          theta,   Nimage);
    174188  gfits_set_bintable_column (&theader, &ftable, "MCAL",           Mcal,    Nimage);
    175189  gfits_set_bintable_column (&theader, &ftable, "SECZ",           secz,    Nimage);
    … …  
    182196  free (crval1);
    183197  free (crval2);
     198  free (theta);
    184199  free (Mcal);
    185200  free (secz);

  • trunk/Ohana/src/libdvo/src/ImageMetadataSelection.c

    r35263 r35416  
    55static off_t Nimage = 0;
    66static Coords mosaic;
     7static Coords field;
    78
    89/* load images based on parameters and region, etc */
    … …  
    1112  image = NULL;
    1213 
     14  /* field defines a frame with 0,0 at the field center, and 1 arcsec / pixel */
     15  field.crpix1 = field.crpix2 = 0.0;
     16  field.cdelt1 = field.cdelt2 = 1.0 / 3600;
     17  field.pc1_1  = field.pc2_2  = 1.0;
     18  field.pc1_2  = field.pc2_1  = 0.0;
     19  field.Npolyterms = 0;
     20  strcpy (field.ctype, "RA---SIN");
     21
    1322  /* mosaic defines a frame with 0,0 at the mosaic center, and 1 arcsec / pixel */
    1423  mosaic.crpix1 = mosaic.crpix2 = 0.0;
    … …  
    4251}
    4352
     53Coords *MatchFieldMetadata (unsigned int imageID) {
     54
     55  int m;
     56
     57  m = match_image_by_ID (image, Nimage, imageID);
     58  if (m == -1) return (NULL);
     59
     60  // if WRP, return the image, otherwise return NULL
     61  // if (strcmp(&image[m].coords.ctype[4], "-WRP")) return NULL;
     62  // return (&image[m].coords);
     63
     64  field.crval1 = image[m].crval1;
     65  field.crval2 = image[m].crval2;
     66  return (&field);
     67}
     68
    4469Coords *MatchMosaicMetadata (unsigned int imageID) {
    4570
    … …  
    4873  m = match_image_by_ID (image, Nimage, imageID);
    4974  if (m == -1) return (NULL);
     75
     76  // if WRP, return the image, otherwise return NULL
     77  // if (strcmp(&image[m].coords.ctype[4], "-WRP")) return NULL;
     78  // return (&image[m].coords);
     79
    5080  mosaic.crval1 = image[m].crval1;
    5181  mosaic.crval2 = image[m].crval2;
     82
     83  mosaic.pc1_1 =  cos(RAD_DEG*image[m].theta);
     84  mosaic.pc1_2 =  sin(RAD_DEG*image[m].theta);
     85  mosaic.pc2_2 =  cos(RAD_DEG*image[m].theta);
     86  mosaic.pc2_1 = -sin(RAD_DEG*image[m].theta);
     87
    5288  return (&mosaic);
    5389}

  • trunk/Ohana/src/libdvo/src/ImageSelection.c

    r35263 r35416  
    9393  int m;
    9494
     95  // mosaic.crval1 = 0;
     96  // mosaic.crval2 = 0;
    9597  m = match_image_subset (image, subset, Nsubset, time, source);
    9698  if (m == -1) return (NULL);
    97   mosaic.crval1 = image[m].coords.crval1;
    98   mosaic.crval2 = image[m].coords.crval2;
    99   return (&mosaic);
     99  // mosaic = image[m].coords.crval1;
     100  // mosaic = image[m].coords.crval2;
     101
     102  // if WRP, return the image, otherwise return NULL
     103  if (strcmp(&image[m].coords.ctype[4], "-WRP")) return NULL;
     104  return (&image[m].coords);
    100105}

  • trunk/Ohana/src/libdvo/src/dbExtractAverages.c

    r34844 r35416  
    141141      value.Flt = average[0].ChiSqPar;
    142142      break;
     143
     144    // XXX case AVE_PM_GROUPS:
     145    // XXX   value.Int = GetProperMotionGroups (average, measure);
     146    // XXX   break;
    143147
    144148    case AVE_TMEAN:
    … …  
    278282} 
    279283
     284// XXX int GetProperMotionGroups (Average *average, Measure *measure) {
     285// XXX   // need the times, excluding ignored detections
     286// XXX   // sort the images
     287// XXX   for (i = 0; i < Ntimes - 1; i++)  {
     288// XXX     if (time[i+1] - time[i] < TRANGE) Ngroup ++;
     289// XXX   }
     290// XXX   return Ngroup;
     291// XXX }
     292// XXX

  • trunk/Ohana/src/libdvo/src/dbExtractMeasures.c

    r35263 r35416  
    1616static int REMOTE_CLIENT = FALSE;
    1717
     18// the following values are calculated together in a single function, eg.,
     19// ApplyTransform() returning Glon & Glat.  for a single measurement, we want to do this
     20// calculation once and save both values in case both are requested (usually both are if
     21// either is)
    1822static int haveGalacticAve = FALSE;
    1923static double GLON_AVE = 0.0;
    … …  
    3135static double ELON_MEAS = 0.0;
    3236static double ELAT_MEAS = 0.0;
     37
     38static int haveMosaicMeas = FALSE;
     39static double XMOS_MEAS = 0.0;
     40static double YMOS_MEAS = 0.0;
     41
     42static int haveFieldMeas = FALSE;
     43static double XFIELD_MEAS = 0.0;
     44static double YFIELD_MEAS = 0.0;
    3345
    3446int dbExtractMeasuresInit (int isRemoteClient) {
    … …  
    6678
    6779int dbExtractMeasuresInitMeas () {
     80  haveMosaicMeas   = FALSE;
    6881  haveGalacticMeas = FALSE;
    6982  haveEclipticMeas = FALSE;
    … …  
    7689  int Nsec;
    7790  dbValue value;
    78   double ra, dec, x, y, dT;
    79 
    80   Coords *mosaic;
     91  double dT;
     92
     93  Coords *mosaic, *fieldc;
     94
    8195  PhotCode *equiv;
    8296
    … …  
    8599
    86100  switch (field->ID) {
    87     case MEAS_MAG: /* magnitudes are already determined above */
    88       equiv = GetPhotcodeEquivbyCode (measure[0].photcode);
    89 
    90       // we return the magnitude for this measure if:
    91       if (field->photcode->type == PHOT_MAG) goto valid_photcode;
     101    case MEAS_MAG: { /* magnitudes are already determined above */
     102      PhotCode *myEquiv = GetPhotcodeEquivbyCode (measure[0].photcode);
     103
     104      // if we request mag:ave, use equiv for photcode
     105      if  (field->photcode->type == PHOT_MAG) {
     106        equiv = myEquiv;
     107        goto valid_photcode;
     108      }
     109
     110      // if we ask for 2MASS_K, etc (REF values), return NAN unless measure->code matches
    92111      if ((field->photcode->type == PHOT_REF) && (measure[0].photcode == field->photcode->code)) goto valid_photcode;
     112
     113      // if we ask for GPC1.g.XY03:rel, etc (DEP values), return NAN unless measure->code matches
    93114      if ((field->photcode->type == PHOT_DEP) && (measure[0].photcode == field->photcode->code)) goto valid_photcode;
    94115
    95       if ((equiv != NULL) && (field->photcode->type == PHOT_SEC) && (equiv[0].code == field->photcode->code)) goto valid_photcode;
     116      // if we ask for g:ave, or other SEC-level values, return the corresponding field
     117      if (field->photcode->type == PHOT_SEC) {
     118        switch (field->magMode) {
     119          // measure-like : return non-NAN if measure.equiv.photcode matches field.photcode
     120          case MAG_INST:
     121          case MAG_CAT:
     122          case MAG_SYS:
     123          case MAG_REL:
     124          case MAG_CAL:
     125          case MAG_APER:
     126          case MAG_APER_INST:
     127          case MAG_KRON:
     128          case MAG_KRON_INST:
     129          case MAG_KRON_ERR:
     130          case MAG_ERR:
     131          case MAG_PHOT_FLAGS:
     132            equiv = myEquiv;
     133            if (equiv && (equiv->code == field->photcode->code)) goto valid_photcode;
     134            break;
     135
     136            // mean-like : return value for the given photcode
     137          case MAG_AVE:
     138          case MAG_REF:
     139          case MAG_CHISQ:
     140          case MAG_AVE_ERR:
     141          case MAG_NCODE:
     142          case MAG_NPHOT:
     143          case MAG_FLUX_PSF:
     144          case MAG_FLUX_PSF_ERR:
     145          case MAG_FLUX_KRON:
     146          case MAG_FLUX_KRON_ERR:
     147            equiv = field->photcode;
     148            goto valid_photcode;
     149            break;
     150          default:
     151            fprintf (stderr, "error");
     152            return value;
     153        }
     154      }
    96155      break;
    97156
    … …  
    182241      }
    183242      break;
     243    }
    184244    case MEAS_RA: /* OK */
    185245      value.Flt = average[0].R - measure[0].dR / 3600.0;
    … …  
    465525      value.Flt = FromShortPixels(measure[0].dRsys);
    466526      break;
    467     case MEAS_XMOSAIC: /* OK */
    468       ra  = average[0].R - measure[0].dR / 3600.0;
    469       dec = average[0].D - measure[0].dD / 3600.0;
    470       mosaic = MatchMosaicMetadata (measure[0].imageID);
    471       if (mosaic == NULL) break;
    472       RD_to_XY (&x, &y, ra, dec, mosaic);
    473       value.Flt = x;
     527
     528    case MEAS_XFIELD: /* offset relative to exposure center in ra,dec space */
     529      if (!haveFieldMeas) {
     530        if (REMOTE_CLIENT) {
     531          fieldc = MatchFieldMetadata (measure[0].imageID);
     532        } else {
     533          fprintf (stderr, "non-parallel Xmos broken\n");
     534          abort();
     535          // fieldc = MatchField (measure[0].t, measure[0].photcode);
     536        }
     537        if (fieldc == NULL) break;
     538        double Rm = average[0].R - measure[0].dR / 3600.0;
     539        double Dm = average[0].D - measure[0].dD / 3600.0;
     540        RD_to_XY (&XFIELD_MEAS, &YFIELD_MEAS, Rm, Dm, fieldc);
     541      }
     542      value.Flt = XFIELD_MEAS;
     543      break;
     544    case MEAS_YFIELD: /* OK */
     545      if (!haveFieldMeas) {
     546        if (REMOTE_CLIENT) {
     547          fieldc = MatchFieldMetadata (measure[0].imageID);
     548        } else {
     549          fprintf (stderr, "non-parallel Xmos broken\n");
     550          abort();
     551          // fieldc = MatchField (measure[0].t, measure[0].photcode);
     552        }
     553        if (fieldc == NULL) break;
     554        double Rm = average[0].R - measure[0].dR / 3600.0;
     555        double Dm = average[0].D - measure[0].dD / 3600.0;
     556        RD_to_XY (&XFIELD_MEAS, &YFIELD_MEAS, Rm, Dm, fieldc);
     557      }
     558      value.Flt = YFIELD_MEAS;
     559      break;
     560
     561    case MEAS_XMOSAIC: /* offset relative to exposure center in camera coords */
     562      if (!haveMosaicMeas) {
     563        if (REMOTE_CLIENT) {
     564          mosaic = MatchMosaicMetadata (measure[0].imageID);
     565        } else {
     566          fprintf (stderr, "non-parallel Xmos broken\n");
     567          abort();
     568          mosaic = MatchMosaic (measure[0].t, measure[0].photcode);
     569        }
     570        if (mosaic == NULL) break;
     571        double Rm = average[0].R - measure[0].dR / 3600.0;
     572        double Dm = average[0].D - measure[0].dD / 3600.0;
     573        RD_to_XY (&XMOS_MEAS, &YMOS_MEAS, Rm, Dm, mosaic);
     574      }
     575      value.Flt = XMOS_MEAS;
    474576      break;
    475577    case MEAS_YMOSAIC: /* OK */
    476       ra  = average[0].R - measure[0].dR / 3600.0;
    477       dec = average[0].D - measure[0].dD / 3600.0;
    478       mosaic = MatchMosaic (measure[0].t, measure[0].photcode);
    479       if (mosaic == NULL) break;
    480       RD_to_XY (&x, &y, ra, dec, mosaic);
    481       value.Flt = y;
     578      if (!haveMosaicMeas) {
     579        if (REMOTE_CLIENT) {
     580          mosaic = MatchMosaicMetadata (measure[0].imageID);
     581        } else {
     582          fprintf (stderr, "non-parallel Xmos broken\n");
     583          abort();
     584          mosaic = MatchMosaic (measure[0].t, measure[0].photcode);
     585        }
     586        if (mosaic == NULL) break;
     587        double Rm = average[0].R - measure[0].dR / 3600.0;
     588        double Dm = average[0].D - measure[0].dD / 3600.0;
     589        RD_to_XY (&XMOS_MEAS, &YMOS_MEAS, Rm, Dm, mosaic);
     590      }
     591      value.Flt = YMOS_MEAS;
    482592      break;
    483593

  • trunk/Ohana/src/libdvo/src/dbFields.c

    r35237 r35416  
    8787    strcpy (code[0].name, "MAG");
    8888    code[0].type = PHOT_MAG;
     89    // the field call 'mag' is only valid for mextract
     90    // it should default to REL, but other types should default to AVE
     91    if (useDefault) {
     92      *mode = MAG_REL;
     93    }
    8994    free (tmpstring);
    9095    return (code);
    … …  
    220225  if (!strcasecmp (fieldName, "YFIX"))           ESCAPE (MEAS_YFIX,           MAG_NONE, OPIHI_FLT);
    221226  if (!strcasecmp (fieldName, "POS_SYS_ERR"))    ESCAPE (MEAS_POS_SYS_ERR,    MAG_NONE, OPIHI_FLT);
     227  if (!strcasecmp (fieldName, "XFIELD"))         ESCAPE (MEAS_XFIELD,         MAG_NONE, OPIHI_FLT);
     228  if (!strcasecmp (fieldName, "YFIELD"))         ESCAPE (MEAS_YFIELD,         MAG_NONE, OPIHI_FLT);
    222229  if (!strcasecmp (fieldName, "XMOSAIC"))        ESCAPE (MEAS_XMOSAIC,        MAG_NONE, OPIHI_FLT);
    223230  if (!strcasecmp (fieldName, "YMOSAIC"))        ESCAPE (MEAS_YMOSAIC,        MAG_NONE, OPIHI_FLT);
    … …  
    254261
    255262  // check for code:mode in photcode name
    256   code = ParsePhotcodeField (fieldName, &mode, MAG_REL);
     263  code = ParsePhotcodeField (fieldName, &mode, MAG_AVE);
    257264  if (code == NULL) {
    258265    gprint (GP_ERR, "unknown field '%s' for measurement table in DVO database\n", fieldName);

  • trunk/Ohana/src/libdvo/src/dvo_catalog.c

    r35102 r35416  
    112112
    113113// init all data, or just catalog data
     114void dvo_averageT_init (AverageTiny *average) {
     115  average->R               = 0;
     116  average->D               = 0;
     117  average->flags           = 0;
     118  average->Nmeasure        = 0;
     119  average->measureOffset   = -1;
     120  average->catID           = 0;
     121}
     122
     123// init all data, or just catalog data
    114124void dvo_secfilt_init (SecFilt *secfilt) {
    115125  secfilt->M           = NAN;
    … …  
    142152 measure->dR        = NAN;
    143153 measure->dD        = NAN;
    144  measure->M         = NAN;
    145154 measure->Mcal      = NAN;
    146155 measure->Map       = NAN;
    … …  
    203212 measure->dbFlags   = 0;
    204213 measure->photFlags = 0;
     214}
     215
     216void dvo_measureT_init (MeasureTiny *measure) {
     217 measure->dR        = NAN;
     218 measure->dD        = NAN;
     219 measure->M         = NAN;
     220 measure->Mcal      = NAN;
     221 measure->dM        = NAN;
     222
     223 measure->airmass   = NAN;
     224 measure->Xccd      = NAN;
     225 measure->Yccd      = NAN;
     226 measure->Xfix      = NAN;
     227 measure->Yfix      = NAN;
     228
     229 measure->t         = 0;
     230 measure->dt        = NAN;
     231 measure->averef    = 0;
     232
     233 measure->imageID   = 0;
     234
     235 measure->dbFlags   = 0;
     236 measure->photFlags = 0;
     237 measure->photcode  = 0;
     238
     239 measure->catID     = 0;
     240
     241 measure->dXccd     = 0;
     242 measure->dYccd     = 0;
     243 measure->dRsys     = 0;
    205244}
    206245
Note: See TracChangeset for help on using the changeset viewer.