IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 35393


Ignore:
Timestamp:
Apr 11, 2013, 2:35:05 PM (13 years ago)
Author:
eugene
Message:

create relphot_images, move relevant code from relphot.c there; add concept of synthetic mags to tie w-band to gri photometry; fix uniphot to handle images with photcodes outside of nominal range

Location:
branches/eam_branches/ipp-20130307/Ohana/src
Files:
2 added
16 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20130307/Ohana/src/addstar/src/MatchHeaders.c

    r35316 r35393  
    6262    if (!strcmp (exttype, "PS1_SV2")) goto keep;
    6363    if (!strcmp (exttype, "PS1_DV3")) {
    64       fprintf (stderr, "exttype, ..%s..\n", exttype);
    6564      goto keep;
    6665    }
    67     fprintf (stderr, "exttype, ..%s..\n", exttype);
    6866    continue;
    6967

  • branches/eam_branches/ipp-20130307/Ohana/src/libdvo/include/dvo.h

    r35316 r35393  
    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
    … …  
    720720
    721721void dvo_average_init (Average *average);
     722void dvo_averageT_init (AverageTiny *average);
    722723void dvo_secfilt_init (SecFilt *secfilt);
    723724void dvo_measure_init (Measure *measure);
     725void dvo_measureT_init (MeasureTiny *measure);
    724726
    725727# endif // DVO_H

  • branches/eam_branches/ipp-20130307/Ohana/src/libdvo/src/dbExtractMeasures.c

    r35342 r35393  
    9292
    9393  Coords *mosaic, *fieldc;
     94
    9495  PhotCode *equiv;
    9596
    … …  
    9899
    99100  switch (field->ID) {
    100     case MEAS_MAG: /* magnitudes are already determined above */
    101       equiv = GetPhotcodeEquivbyCode (measure[0].photcode);
    102 
    103       // we return the magnitude for this measure if:
    104       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
    105111      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
    106114      if ((field->photcode->type == PHOT_DEP) && (measure[0].photcode == field->photcode->code)) goto valid_photcode;
    107115
    108       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      }
    109155      break;
    110156
    … …  
    195241      }
    196242      break;
     243    }
    197244    case MEAS_RA: /* OK */
    198245      value.Flt = average[0].R - measure[0].dR / 3600.0;

  • branches/eam_branches/ipp-20130307/Ohana/src/libdvo/src/dbFields.c

    r35342 r35393  
    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);
    … …  
    256261
    257262  // check for code:mode in photcode name
    258   code = ParsePhotcodeField (fieldName, &mode, MAG_REL);
     263  code = ParsePhotcodeField (fieldName, &mode, MAG_AVE);
    259264  if (code == NULL) {
    260265    gprint (GP_ERR, "unknown field '%s' for measurement table in DVO database\n", fieldName);

  • branches/eam_branches/ipp-20130307/Ohana/src/libdvo/src/dvo_catalog.c

    r35102 r35393  
    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

  • branches/eam_branches/ipp-20130307/Ohana/src/relphot/Makefile

    r34844 r35393  
    4040$(SRC)/load_catalogs.$(ARCH).o   \
    4141$(SRC)/load_images.$(ARCH).o     \
     42$(SRC)/synthetic_mags.$(ARCH).o  \
    4243$(SRC)/plot_scatter.$(ARCH).o    \
    4344$(SRC)/plotstuff.$(ARCH).o       \
    4445$(SRC)/reload_catalogs.$(ARCH).o \
    4546$(SRC)/relphot.$(ARCH).o         \
     47$(SRC)/relphot_images.$(ARCH).o  \
    4648$(SRC)/relphot_objects.$(ARCH).o         \
    4749$(SRC)/select_images.$(ARCH).o   \
    … …  
    6971$(SRC)/args.$(ARCH).o            \
    7072$(SRC)/help.$(ARCH).o            \
     73$(SRC)/synthetic_mags.$(ARCH).o  \
    7174$(SRC)/plotstuff.$(ARCH).o       \
    7275$(SRC)/liststats.$(ARCH).o       \

  • branches/eam_branches/ipp-20130307/Ohana/src/relphot/include/relphot.h

    r35104 r35393  
    44# include <signal.h>
    55# include <pthread.h>
     6
     7# define MARKTIME(MSG,...) {                    \
     8    gettimeofday (&stopTimer, (void *) NULL);   \
     9    float dtime = DTIME (stopTimer, startTimer);        \
     10    fprintf (stderr, MSG, __VA_ARGS__); }
     11
     12# define INITTIME \
     13  struct timeval startTimer, stopTimer; \
     14  gettimeofday (&startTimer, (void *) NULL);
    615
    716/* # define GRID_V1 */
    … …  
    193202int    UpdateAverages;
    194203int    ApplyOffsets;
     204int    SyntheticPhotometry;
    195205
    196206char  *PhotcodeList;
    … …  
    353363void          write_coords        PROTO((Header *header, Coords *coords));
    354364int relphot_objects (int hostID, char *hostpath);
     365int relphot_images (void);
    355366
    356367void relphot_usage (void);
    … …  
    390401
    391402void SetMrelInfoInit (SetMrelInfo *results, int allocLists);
     403
     404int init_synthetic_mags ();
     405int add_synthetic_mags (AverageTiny *average, SecFilt *secfilt, MeasureTiny *measure, off_t *Nmeasure, off_t *Nm);

  • branches/eam_branches/ipp-20130307/Ohana/src/relphot/src/MosaicOps.c

    r34088 r35393  
    152152}
    153153
    154 # define MARKTIME(MSG,...) { \
    155   float dtime; \
    156   gettimeofday (&stopTimer, (void *) NULL); \
    157   dtime = DTIME (stopTimer, startTimer); \
    158   fprintf (stderr, MSG, __VA_ARGS__); }
    159 
    160154/* find mosaic frames (unique time periods) (NOTE : we do NOT require matching photcodes...)
    161155   this function will also identify the images NOT in the subset which belong to a selected mosaic
    … …  
    167161  char *pname;
    168162
    169   struct timeval startTimer, stopTimer;
    170 
    171163  if (!MOSAIC_ZEROPT) return;
    172164
    173   gettimeofday (&startTimer, (void *) NULL);
     165  INITTIME;
    174166
    175167  /* a 'mosaic' in relphot is (unlike relastro) a virtual concept: there is no
    … …  
    264256    }
    265257  }
    266   fprintf (stderr, "%d total images, %d overlap skyregion, %d do not overlap, but match overlapping mosaics\n", (int) Nimage, (int) Nskip, (int) Nmark);
     258  fprintf (stderr, "%d total images, %d overlap skyregion & match selection criteria, %d do not overlap, but match overlapping mosaics\n", (int) Nimage, (int) Nskip, (int) Nmark);
    267259  MARKTIME("find unselected images matching selected mosaics: %f sec\n", dtime);
    268260

  • branches/eam_branches/ipp-20130307/Ohana/src/relphot/src/args.c

    r35104 r35393  
    251251    remove_argument (N, &argc, argv);
    252252    MOSAIC_ZEROPT = TRUE;
     253    if (!strcasecmp (MOSAICNAME, "none")) {
     254      fprintf (stderr, "mosaic astrometry selected by MOSAICNAME not defined\n");
     255      exit (2);
     256    }
    253257  }
    254258
    … …  
    315319    remove_argument (N, &argc, argv);
    316320    DophotSelect = TRUE;
     321  }
     322
     323  SyntheticPhotometry = FALSE;
     324  if ((N = get_argument (argc, argv, "-synthphot"))) {
     325    remove_argument (N, &argc, argv);
     326    SyntheticPhotometry = TRUE;
     327    init_synthetic_mags();
    317328  }
    318329
    … …  
    516527  }
    517528
     529  SyntheticPhotometry = FALSE;
     530  if ((N = get_argument (argc, argv, "-synthphot"))) {
     531    remove_argument (N, &argc, argv);
     532    SyntheticPhotometry = FALSE;
     533    init_synthetic_mags();
     534  }
     535
    518536  ImagSelect = FALSE;
    519537  if ((N = get_argument (argc, argv, "-instmag"))) {

  • branches/eam_branches/ipp-20130307/Ohana/src/relphot/src/bcatalog.c

    r34405 r35393  
    179179        REALLOCATE (subcatalog[0].measureT, MeasureTiny, NMEASURE);
    180180      }
     181      myAssert (Nmeasure < NMEASURE, "realloc failure");
    181182    }
    182183
    … …  
    204205      continue;
    205206    }
     207
     208    // for w-band photometry (& other cases?) convert gri photometry
     209    // to a synthetic w-band magnitude
     210    add_synthetic_mags (&subcatalog[0].averageT[Naverage],
     211                        &subcatalog[0].secfilt[Nsecfilt*Naverage],
     212                        &subcatalog[0].measureT[Nmeasure],
     213                        &Nmeasure,
     214                        &Nm);
     215    if (Nmeasure == NMEASURE) {
     216      NMEASURE += 1000;
     217      REALLOCATE (subcatalog[0].measureT, MeasureTiny, NMEASURE);
     218    }
     219    myAssert (Nmeasure < NMEASURE, "realloc failure");
    206220
    207221    subcatalog[0].averageT[Naverage].Nmeasure = Nm;

  • branches/eam_branches/ipp-20130307/Ohana/src/relphot/src/global_stats.c

    r33963 r35393  
    11# include "relphot.h"
    2 
    3 # define MARKTIME(MSG,...) { \
    4   float dtime; \
    5   gettimeofday (&stopTimer, (void *) NULL); \
    6   dtime = DTIME (stopTimer, startTimer); \
    7   fprintf (stderr, MSG, __VA_ARGS__); }
    82
    93void global_stats (Catalog *catalog, int Ncatalog, FlatCorrectionTable *flatcorr) {
    … …  
    115  StatType stN, stX, stS, imN, imX, imM, imD, msM, msX, msN, msD;
    126
    13   // struct timeval startTimer, stopTimer;
     7  // INITTIME;
    148
    159  fprintf (stderr, "\n");
    1610  fprintf (stderr, "STATS            median     mean    sigma      min      max   Nmeas\n");
    17 
    18   // gettimeofday (&startTimer, (void *) NULL);
    1911
    2012  int Ns;

  • branches/eam_branches/ipp-20130307/Ohana/src/relphot/src/load_catalogs.c

    r34642 r35393  
    148148
    149149    char tmpline[1024];
    150     if (VERBOSE)       { snprintf (tmpline, 1024, "%s -v",              command);                    strcpy (command, tmpline); }
    151     if (VERBOSE2)      { snprintf (tmpline, 1024, "%s -vv",             command);                    strcpy (command, tmpline); }
    152     if (RESET)         { snprintf (tmpline, 1024, "%s -reset",          command);                    strcpy (command, tmpline); }
    153     if (RESET_ZEROPTS) { snprintf (tmpline, 1024, "%s -reset-zpts",     command);                    strcpy (command, tmpline); }
    154     if (!KEEP_UBERCAL) { snprintf (tmpline, 1024, "%s -reset-ubercal",  command);                    strcpy (command, tmpline); }
    155     if (DophotSelect)  { snprintf (tmpline, 1024, "%s -dophot %d",      command, DophotValue);       strcpy (command, tmpline); }
    156     if (ImagSelect)    { snprintf (tmpline, 1024, "%s -instmag %f %f",  command, ImagMin, ImagMax);  strcpy (command, tmpline); }
    157     if (MaxDensityUse) { snprintf (tmpline, 1024, "%s -max-density %f", command, MaxDensityValue);   strcpy (command, tmpline); }
     150    if (VERBOSE)             { snprintf (tmpline, 1024, "%s -v",              command);                   strcpy (command, tmpline); }
     151    if (VERBOSE2)            { snprintf (tmpline, 1024, "%s -vv",             command);                   strcpy (command, tmpline); }
     152    if (RESET)               { snprintf (tmpline, 1024, "%s -reset",          command);                   strcpy (command, tmpline); }
     153    if (RESET_ZEROPTS)       { snprintf (tmpline, 1024, "%s -reset-zpts",     command);                   strcpy (command, tmpline); }
     154    if (!KEEP_UBERCAL)       { snprintf (tmpline, 1024, "%s -reset-ubercal",  command);                   strcpy (command, tmpline); }
     155    if (DophotSelect)        { snprintf (tmpline, 1024, "%s -dophot %d",      command, DophotValue);      strcpy (command, tmpline); }
     156    if (ImagSelect)          { snprintf (tmpline, 1024, "%s -instmag %f %f",  command, ImagMin, ImagMax); strcpy (command, tmpline); }
     157    if (MaxDensityUse)       { snprintf (tmpline, 1024, "%s -max-density %f", command, MaxDensityValue);  strcpy (command, tmpline); }
     158    if (SyntheticPhotometry) { snprintf (tmpline, 1024, "%s -synthphot",      command);                   strcpy (command, tmpline); }
     159
    158160    if (TimeSelect) {
    159161      char *tstart = ohana_sec_to_date (TSTART);

  • branches/eam_branches/ipp-20130307/Ohana/src/relphot/src/load_images.c

    r33963 r35393  
    11# include "relphot.h"
    2 
    3 # define MARKTIME(MSG,...) { \
    4   float dtime; \
    5   gettimeofday (&stop, (void *) NULL); \
    6   dtime = DTIME (stop, start); \
    7   fprintf (stderr, MSG, __VA_ARGS__); }
    82
    93// This function generates a subset of the images based on selections.  Input db has already
    … …  
    148  off_t      Nimage, Nsubset, Nchar;
    159  off_t     *LineNumber;
    16   struct timeval start, stop;
    1710  char *inSubset;
    1811
    … …  
    2013  SkyList *skylist = NULL;
    2114
    22   gettimeofday (&start, (void *) NULL);
     15  INITTIME;
    2316
    2417  // load the current sky table (layout of all SkyRegions)
    … …  
    9386  }
    9487
    95 # if (1)
     88# if (0)
    9689  FILE *ftest = fopen ("skydump.dat", "w");
    9790  int i;

  • branches/eam_branches/ipp-20130307/Ohana/src/relphot/src/relphot.c

    r33963 r35393  
    11# include "relphot.h"
    22
    3 # define MARKTIME(MSG,...) {                    \
    4     float dtime;                                \
    5     gettimeofday (&stop, (void *) NULL);        \
    6     dtime = DTIME (stop, start);                \
    7     fprintf (stderr, MSG, __VA_ARGS__); }
    8 
    93int main (int argc, char **argv) {
    10 
    11   int i, status, Ncatalog;
    12   Catalog *catalog = NULL;
    13   FITS_DB db;
    14   struct timeval start, stop;
    15   SkyList *skylist = NULL;
    16 
    17   gettimeofday (&start, (void *) NULL);
    184
    195  /* get configuration info, args */
    … …  
    3622  }
    3723
    38   /* register database handle with shutdown procedure */
    39   set_db (&db);
    40   db.mode   = dvo_catalog_catmode (CATMODE);
    41   db.format = dvo_catalog_catformat (CATFORMAT);
    42 
    43   /* lock and load the image db table */
    44   status = dvo_image_lock (&db, ImageCat, 60.0, (UPDATE ? LCK_XCLD : LCK_SOFT));
    45   if (!status && UPDATE) {
    46     fprintf (stderr, "error\n");
    47     Shutdown ("ERROR: failure to lock image catalog %s", db.filename);
    48   }
    49   // if the file is missing, db.dbstate will have a value of either:
    50   // LCK_EMPTY (if UPDATE) or LCK_MISSING (if !UPDATE)
    51   if ((db.dbstate == LCK_EMPTY) || (db.dbstate == LCK_MISSING)) {
    52 
    53     // 25.0 is the nominal zero point for measurements in the database
    54     // for all measurements except PHOT_REF:
    55     // measure.M contains: -2.5*log(DN) + 2.5*log(exptime) + 25.0
    56     // for measurements using PHOT_REF:
    57     // measure.M contains: -2.5*log(DN) + 2.5*log(exptime) + True Zero Point (but exptime and ZP are known)
    58     dvo_image_create (&db, 25.0);
    59 
    60     // XXX why are we running relphot on an empty database?
    61   } else {
    62     if (!dvo_image_load (&db, VERBOSE, FALSE)) Shutdown ("can't read image catalog %s", db.filename);
    63   }
    64   MARKTIME("-- load image data: %f sec\n", dtime);
    65 
    66   /* load regions and images based on specified sky patch and/or catalog */
    67   skylist = load_images (&db, UserCatalog, &UserPatch);
    68   MARKTIME("-- load images: %f sec\n", dtime);
    69 
    70   /* unlock, if we can (else, unlocked below) */
    71   if (!UPDATE) dvo_image_unlock (&db);
    72 
    73   // load the flat correction table (if defined)
    74   char flatcorrfile[256];
    75   sprintf (flatcorrfile, "%s/flatcorr.fits", CATDIR);
    76   FlatCorrectionTable *flatcorr = FlatCorrectionLoad (flatcorrfile, VERBOSE);
    77 
    78   if (NLOOP > 0) {
    79     /* load catalog data from region files (hostID is 0 since we are not a client */
    80     catalog = load_catalogs (skylist, &Ncatalog, 0, NULL);
    81     MARKTIME("-- load catalog data: %f sec\n", dtime);
    82  
    83     /* add in a loop over the catalogs calling dvo_catalog_chipcoords */
    84 
    85     /* match measurements with images, mosaics */
    86     initImageBins  (catalog, Ncatalog, TRUE);
    87     MARKTIME("-- make image bins: %f sec\n", dtime);
    88 
    89     initMosaicBins (catalog, Ncatalog, TRUE);
    90     initGridBins   (catalog, Ncatalog);
    91     initMrel (catalog, Ncatalog);
    92 
    93     findImages (catalog, Ncatalog, TRUE);
    94     MARKTIME("-- set up image indexes: %f sec\n", dtime);
    95 
    96     findMosaics (catalog, Ncatalog, TRUE);  /* also sets Grid values */
    97     MARKTIME("-- set up mosaic indexes: %f sec\n", dtime);
    98 
    99     SAVEPLOT = FALSE;
    100 
    101     setExclusions (catalog, Ncatalog, TRUE);
    102 
    103     global_stats (catalog, Ncatalog, flatcorr);
    104 
    105     if (PLOTSTUFF) {
    106       plot_star_coords (catalog, Ncatalog);
    107       // plot_mosaic_fields (catalog);
    108     }
    109 
    110     // if we are measuring the flat-field correction grid, we need to perform a number of iterations first:
    111     if (USE_GRID) {
    112       int star_toofew;
    113 
    114       // until we finish the grid analysis, do not reject stars out-of-hand based on ID_STAR_FEW
    115       // XXX this is kind of poor: need to have a better distinctions about STAR_BAD in setMrel vs getMrel
    116       star_toofew = STAR_TOOFEW;
    117       STAR_TOOFEW = 0;
    118       for (i = 0; i < NGRID; i++) {
    119         STAR_BAD = ID_STAR_POOR;
    120         setMrel  (catalog, Ncatalog, flatcorr);
    121         STAR_BAD  = ID_STAR_POOR | ID_STAR_FEW;
    122         setMgrid (catalog, flatcorr);
    123       }
    124       STAR_BAD  = ID_STAR_POOR | ID_STAR_FEW;
    125       STAR_TOOFEW = star_toofew;
    126     }
    127 
    128     /* determine fit values */
    129     for (i = 0; i < NLOOP; i++) {
    130       if (PLOTSTUFF) {
    131         plot_scatter (catalog, Ncatalog, flatcorr);
    132       }
    133       setMrel  (catalog, Ncatalog, flatcorr); // threaded
    134       if (PLOTSTUFF) {
    135         plot_scatter (catalog, Ncatalog, flatcorr);
    136       }
    137       setMcal  (catalog, FALSE, flatcorr);
    138       setMmos  (catalog, FALSE, flatcorr);
    139       setMgrid (catalog, flatcorr);
    140       MARKTIME("-- set Mrel, Mcal, Mmos, Mgrid : %f sec\n", dtime);
    141    
    142       if (PLOTSTUFF) {
    143         plot_scatter (catalog, Ncatalog, flatcorr);
    144         plot_grid (catalog, flatcorr);
    145         plot_mosaics ();
    146         plot_images ();
    147         plot_stars (catalog, Ncatalog);
    148         plot_chisq (catalog, Ncatalog);
    149       }
    150       // if (i < NLOOP - 1) rationalize_mosaics (catalog, Ncatalog);
    151       // if (i % 6 == 1) rationalize_images ();
    152 
    153       // NOTE : in the past, I was not iterating enough before cleaning.  make sure we do
    154       // at least 8 loops first -- that should get the systematic errors down to the ~1%
    155       // level, even in cases where we have an even split between photometric data and
    156       // data with 1 mag of extinction.
    157 
    158       if ((i > 8) && (i % 8 == 2)) clean_measures (catalog, Ncatalog, FALSE, flatcorr);
    159       if ((i > 8) && (i % 8 == 3)) clean_stars (catalog, Ncatalog);
    160       if ((i > 8) && (i % 8 == 5)) clean_mosaics ();
    161       if ((i > 8) && (i % 8 == 5)) clean_images ();
    162 
    163       // if ((i == 1) || (i == 5) || (i ==  9) || (i == 13)) clean_measures (catalog, Ncatalog, FALSE);
    164       // if ((i == 2) || (i == 6) || (i == 10) || (i == 14)) clean_stars (catalog, Ncatalog);
    165       // if ((i == 4) || (i == 8) || (i == 12) || (i == 16)) clean_mosaics ();
    166       // if ((i == 4) || (i == 8) || (i == 12) || (i == 16)) clean_images ();
    167       if (i % 3 == 2) global_stats (catalog, Ncatalog, flatcorr);
    168       MARKTIME("-- finished loop %d: %f sec\n", i, dtime);
    169     }
    170 
    171     if (PLOTSTUFF) {
    172       plot_scatter (catalog, Ncatalog, flatcorr);
    173       plot_grid (catalog, flatcorr);
    174       plot_mosaics ();
    175       plot_images ();
    176       plot_stars (catalog, Ncatalog);
    177       plot_chisq (catalog, Ncatalog);
    178     }
    179  
    180     if (USE_GRID) dump_grid ();
    181 
    182     /* set Mcal & Mmos for bad images */
    183     setMcal  (catalog, TRUE, flatcorr);
    184     setMmos  (catalog, TRUE, flatcorr);
    185     MARKTIME("-- finalize Mcal values: %f sec\n", dtime);
    186 
    187     setMcalFinal (); // copy per-mosaic calibrations to the images
    188 
    189     if (SAVE_IMAGE_UPDATES) {
    190 
    191       FITS_DB dbX;
    192       off_t *LineNumber;
    193       Image *image, *subset;
    194       off_t Nimage, Nsubset, i, Nx;
    195       char filename[1024];
    196 
    197       gfits_db_init (&dbX);
    198       dbX.lockstate = LCK_XCLD;
    199       dbX.timeout   = 60.0;
    200       dbX.mode      = db.mode;
    201       dbX.format    = db.format;
    202 
    203       snprintf (filename, 1024, "%s.bck", ImageCat);
    204       if (!gfits_db_lock (&dbX, filename)) {
    205         fprintf (stderr, "can't lock backup image image catalog\n");
    206         return (FALSE);
    207       }
    208    
    209       // apply the changes from the image subset to the full table:
    210 
    211       // convert database table to internal structure (binary to Image)
    212       // 'image' points to the same memory as db->ftable->buffer
    213       image = gfits_table_get_Image (&db.ftable, &Nimage, &db.swapped);
    214       if (!image) {
    215         fprintf (stderr, "ERROR: failed to read images\n");
    216         exit (2);
    217       }
    218       gfits_scan (db.ftable.header, "NAXIS1", OFF_T_FMT, 1,  &Nx);
    219 
    220       subset = getimages (&Nsubset, &LineNumber);
    221       for (i = 0; i < Nsubset; i++) {
    222         if (LineNumber[i] == -1) continue;
    223         memcpy (&image[LineNumber[i]], &subset[i], Nx);
    224       }
    225 
    226       // copy Images.dat data (all or only the subset / vtable elements?)  I think I need to dump
    227       // the entire Image table, but with the updates in place I think this says: re-work the
    228       // ftable/vtable usage in this program to be more sensible...
    229       gfits_copy_header (&db.header,  &dbX.header);
    230       gfits_copy_matrix (&db.matrix,  &dbX.matrix);
    231       gfits_copy_header (&db.theader, &dbX.theader);
    232       gfits_copy_ftable (&db.ftable,  &dbX.ftable);
    233 
    234       dbX.ftable.header = &dbX.theader;
    235       dbX.virtual = FALSE;
    236 
    237       // save Images.dat using the copied structure
    238       if (UPDATE_CATFORMAT) {
    239         // ensure the db format is updated
    240         dbX.format = dvo_catalog_catformat (UPDATE_CATFORMAT);
    241         char photcodeFile[1024];
    242         sprintf (photcodeFile, "%s/Photcodes.dat", CATDIR);
    243         SavePhotcodesFITS (photcodeFile);
    244       }
    245       dvo_image_save (&dbX, VERBOSE);
    246       dvo_image_unlock (&dbX);
    247       MARKTIME("-- save Image.dat.bck: %f sec\n", dtime);
    248     }
    249 
    250     /* at this point, we have correct cal coeffs in the image/mosaic structures */
    251     for (i = 0; i < Ncatalog; i++) {
    252       free_tiny_values (&catalog[i]);
    253       dvo_catalog_free (&catalog[i]);
    254     }
    255     freeImageBins (Ncatalog, TRUE);
    256     freeMosaicBins (Ncatalog, TRUE);
    257     freeGridBins (Ncatalog);
    258 
    259   } // (NLOOP > 0) : this loop determines the offsets per chip
    260 
    261   reload_images (&db);
    262 
    263   // only change the real database files if -update is requested
    264   if (!UPDATE) exit (0);
    265  
    266   /* load catalog data from region files, update Mrel include all data */
    267   reload_catalogs (skylist, flatcorr, 0, NULL);
    268   MARKTIME("-- updated all catalogs: %f sec\n", dtime);
    269 
    270   if (UPDATE_CATFORMAT) {
    271     // ensure the db format is updated
    272     db.format = dvo_catalog_catformat (UPDATE_CATFORMAT);
    273   }
    274   dvo_image_update (&db, VERBOSE);
    275   dvo_image_unlock (&db);
     24  relphot_images ();
    27625
    27726  exit (0);

  • branches/eam_branches/ipp-20130307/Ohana/src/relphot/src/select_images.c

    r34642 r35393  
    1717off_t getRegionStartByRA (double R, double *Rref, off_t Nregions);
    1818
    19 # define MARKTIME(MSG,...) { \
    20   float dtime; \
    21   gettimeofday (&stop, (void *) NULL); \
    22   dtime = DTIME (stop, start); \
    23   fprintf (stderr, MSG, __VA_ARGS__); }
    24 
    2519Image *select_images (SkyList *skylist, Image *timage, off_t Ntimage, char *inSubset, off_t **LineNumber, off_t *Nimage, SkyRegion *region) {
    2620 
    … …  
    3125  double Ri[5], Di[5], Xi[5], Yi[5];
    3226  Coords tcoords;
    33   struct timeval start, stop;
    3427 
    3528  if (skylist[0].Nregions < 1) {
    … …  
    4033  }
    4134
    42   gettimeofday (&start, (void *) NULL);
     35  INITTIME;
    4336
    4437  // the comparison is made in the catalog local projection. below we set crval1,2
    … …  
    8275  }
    8376
    84 # if (1)
     77# if (0)
    8578  // XXX quick test of the dec bands:
    8679  FILE *fout = fopen ("dec.bands.dat", "w");
    … …  
    132125    // this adds 1.3 sec for 3M images
    133126    if (!FindMosaicForImage (timage, Ntimage, i)) {
    134       fprintf (stderr, "cannot find mosaic for "OFF_T_FMT"\n", i);
     127      if (VERBOSE2) fprintf (stderr, "cannot find mosaic for "OFF_T_FMT"\n", i);
    135128      continue;
    136129    }

  • branches/eam_branches/ipp-20130307/Ohana/src/uniphot/src/match_flatcorr_to_images.c

    r33654 r35393  
    4343  for (i = 0; i < Nimage; i++) {
    4444    if (!image[i].photcode) continue; // skip PHU images
     45    if (image[i].photcode > maxCode) continue; // skip PHU images
    4546
    4647    int found = FALSE;
Note: See TracChangeset for help on using the changeset viewer.