IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 39272 for trunk


Ignore:
Timestamp:
Dec 14, 2015, 6:14:36 PM (11 years ago)
Author:
eugene
Message:

setphot is nearly working

Location:
trunk/Ohana/src/uniphot
Files:
9 edited

Legend:

Unmodified
Added
Removed
  • trunk/Ohana/src/uniphot/Makefile

    r39262 r39272  
    3737
    3838# install: $(DESTBIN)/uniphot $(DESTBIN)/setfwhm $(DESTBIN)/setphot $(DESTBIN)/setphot_client $(DESTBIN)/setposangle $(DESTBIN)/setposangle_client $(DESTBIN)/setastrom $(DESTBIN)/setastrom_client $(DESTBIN)/setgalmodel $(DESTBIN)/setgalmodel_client $(DESTBIN)/fiximids $(DESTBIN)/fiximids_client $(DESTBIN)/fixstkids $(DESTBIN)/fixstkids_client $(DESTBIN)/ckids $(DESTBIN)/ckids_client
    39 install: $(DESTBIN)/uniphot $(DESTBIN)/setfwhm $(DESTBIN)/setposangle $(DESTBIN)/setposangle_client $(DESTBIN)/setastrom $(DESTBIN)/setastrom_client $(DESTBIN)/setgalmodel $(DESTBIN)/setgalmodel_client $(DESTBIN)/fiximids $(DESTBIN)/fiximids_client $(DESTBIN)/fixstkids $(DESTBIN)/fixstkids_client $(DESTBIN)/ckids $(DESTBIN)/ckids_client
     39install: $(DESTBIN)/uniphot $(DESTBIN)/setfwhm $(DESTBIN)/setphot $(DESTBIN)/setposangle $(DESTBIN)/setposangle_client $(DESTBIN)/setastrom $(DESTBIN)/setastrom_client $(DESTBIN)/setgalmodel $(DESTBIN)/setgalmodel_client $(DESTBIN)/fiximids $(DESTBIN)/fiximids_client $(DESTBIN)/fixstkids $(DESTBIN)/fixstkids_client $(DESTBIN)/ckids $(DESTBIN)/ckids_client
    4040
    4141UNIPHOT = \
     
    6868$(SRC)/match_zpts_to_images.$(ARCH).o       \
    6969$(SRC)/match_flatcorr_to_images.$(ARCH).o           \
     70$(SRC)/match_camcorr_to_images.$(ARCH).o            \
     71$(SRC)/cam_zpt_correction.$(ARCH).o         \
    7072$(SRC)/update_dvo_setphot.$(ARCH).o \
    7173$(SRC)/update_catalog_setphot.$(ARCH).o \
  • trunk/Ohana/src/uniphot/include/setphot.h

    r39262 r39272  
    3535  int NyCCD;    // number of pixels
    3636
     37  e_time *tstart;
     38  e_time *tstop;
     39
     40  Header *phu;
    3741  Matrix **matrix; // allocate an array of pointers
     42  Header **header; // allocate an array of pointers
     43 
    3844  // index = ix + iy*Nx + filter*Nchips + dir*Ngroup
    3945} CamCorrection;
     
    4753char        *IMAGES;
    4854char        *SINGLE_CPT;
     55char        *CAM_FILE;
    4956int          VERBOSE;
    5057int          RESET;
     
    8996Image        *load_images_setphot    PROTO((FITS_DB *db, off_t *Nimage));
    9097int           match_zpts_to_images   PROTO((Image *image, off_t Nimage, ZptTable *zpts, int Nzpts));
    91 int           update_dvo_setphot     PROTO((Image *image, off_t Nimage, FlatCorrectionTable *flatcorr));
    92 void          update_catalog_setphot PROTO((Catalog *catalog, Image *image, off_t *index, off_t Nimage, FlatCorrectionTable *flatcorr));
     98int           update_dvo_setphot     PROTO((Image *image, off_t Nimage, CamCorrection *camcorr));
     99void          update_catalog_setphot PROTO((Catalog *catalog, Image *image, off_t *index, off_t Nimage, CamCorrection *camcorr));
    93100void          update_catalog_setphot_client PROTO((Catalog *catalog, ImageSubset *image, off_t *index, off_t Nimage, FlatCorrectionTable *flatcorr));
    94101
     
    100107float         apply_zpt_offset            PROTO((short code));
    101108
    102 int repair_catalog_by_objID (Catalog *catalog);
     109int           repair_catalog_by_objID     PROTO((Catalog *catalog));
     110
     111void           free_cam_correction          PROTO((CamCorrection *cam));
     112CamCorrection *alloc_cam_correction         PROTO((int Nx, int Ny, int Nfilter, int Nseason));
     113CamCorrection *load_cam_correction          PROTO((char *filename));
     114int            save_cam_correction          PROTO((CamCorrection *cam, char *filename));
     115// int            get_cam_correction        PROTO((int chipID, int filter, float Xccd, float Yccd, double *dX, double *dY));
     116float          get_cam_correction           PROTO((CamCorrection *cam, int flat_id, float Xccd, float Yccd));
     117CamCorrection *merge_flatcorr_and_camcorr   PROTO((FlatCorrectionTable *flatcorr, CamCorrection *camcorr));
     118
     119int            match_camcorr_to_images      PROTO((Image *image, off_t Nimage, CamCorrection *camcorr));
  • trunk/Ohana/src/uniphot/src/cam_zpt_correction.c

    r39267 r39272  
    11# include "setphot.h"
     2
     3// hard-coded values
     4char filters[5][3] = {"g", "r", "i", "z", "y"};
    25
    36// camera systematic (photoflat) correction functions:
     
    69// note that we are allocating pointers for XY00 - XY77, but only the octal elements are set (and not the 00,07,70,77 ones)
    710
    8 static CamCorrection *cam = NULL;
    9 
    10 CamCorrection *get_cam_correction_ptr () {
    11   return cam;
    12 }
    13 
    1411void free_cam_correction (CamCorrection *cam) {
    1512
    1613  if (!cam) return;
    1714
     15  int i;
    1816  for (i = 0; i < cam->Nvalues; i++) {
    1917    FREE (cam->matrix[i]);
    20   }
     18    FREE (cam->header[i]);
     19  }
     20  FREE (cam->phu);
     21
     22  free (cam->tstart);
     23  free (cam->tstop);
    2124
    2225  free (cam->matrix);
     26  free (cam->header);
    2327  free (cam);
    2428}
     
    3034  ALLOCATE (cam, CamCorrection, 1);
    3135
     36  cam->phu = NULL;
     37  ALLOCATE (cam->header, Header *, Nseason*Nfilter*Nx*Ny);
    3238  ALLOCATE (cam->matrix, Matrix *, Nseason*Nfilter*Nx*Ny);
     39  ALLOCATE (cam->tstart, e_time,   Nseason);
     40  ALLOCATE (cam->tstop,  e_time,   Nseason);
    3341
    3442  cam->Nx      = Nx;       // for gpc1, should be 8
     
    4250
    4351  ALLOCATE (cam->matrix, Matrix *, cam->Nvalues);
     52
     53  int i;
    4454  for (i = 0; i < cam->Nvalues; i++) {
    4555    cam->matrix[i] = NULL;
     
    4959}
    5060
    51 int load_cam_correction (char *filename) {
    52 
    53   int i, ix, iy, dir, filter;
    54   Header header;
     61CamCorrection *load_cam_correction (char *filename) {
     62
     63  int i, ix, iy, season, filter;
     64
     65  Header *phu = gfits_alloc_header();
    5566
    5667  /* open file for input */
     
    6273
    6374  // read the PHU header
    64   if (!gfits_load_header (f, &header)) return FALSE;
     75  if (!gfits_load_header (f, phu)) return FALSE;
    6576 
    66   int Nbytes = gfits_data_size (&header);
     77  int Nbytes = gfits_data_size (phu);
    6778  fseeko (f, Nbytes, SEEK_CUR);
    6879
    69   int Ndir, Nfilter, Nx, Ny, dX, dY, NxCCD, NyCCD;
    70   if (!gfits_scan (&header, "NSEASON",  "%d", 1, &Nseason)) return FALSE;
    71   if (!gfits_scan (&header, "NFILTER",  "%d", 1, &Nfilter)) return FALSE;
    72   if (!gfits_scan (&header, "NX",       "%d", 1, &Nx))      return FALSE;
    73   if (!gfits_scan (&header, "NY",       "%d", 1, &Ny))      return FALSE;
    74   if (!gfits_scan (&header, "DX",       "%d", 1, &dX))      return FALSE;
    75   if (!gfits_scan (&header, "DY",       "%d", 1, &dY))      return FALSE;
    76 
    77   if (!gfits_scan (&header, "NX_CHIP", "%d", 1, &NxCCD))   return FALSE;
    78   if (!gfits_scan (&header, "NY_CHIP", "%d", 1, &NyCCD))   return FALSE;
     80  int Nseason, Nfilter, Nx, Ny, dX, dY, NxCCD, NyCCD;
     81  if (!gfits_scan (phu, "NSEASON", "%d", 1, &Nseason)) return FALSE;
     82  if (!gfits_scan (phu, "NFILTER", "%d", 1, &Nfilter)) return FALSE;
     83  if (!gfits_scan (phu, "NX",      "%d", 1, &Nx))      return FALSE;
     84  if (!gfits_scan (phu, "NY",      "%d", 1, &Ny))      return FALSE;
     85  if (!gfits_scan (phu, "DX",      "%d", 1, &dX))      return FALSE;
     86  if (!gfits_scan (phu, "DY",      "%d", 1, &dY))      return FALSE;
     87  if (!gfits_scan (phu, "NX_CHIP", "%d", 1, &NxCCD))   return FALSE;
     88  if (!gfits_scan (phu, "NY_CHIP", "%d", 1, &NyCCD))   return FALSE;
    7989
    8090  CamCorrection *cam = alloc_cam_correction (Nx, Ny, Nfilter, Nseason);
     91  cam->phu = phu;
    8192
    8293  cam->dX      = dX;       // superpixel sampling, x-dir
     
    8697  cam->NyCCD   = NyCCD;    // pixels per chip
    8798
     99  for (i = 0; i < Nseason; i++) {
     100    char field[256];
     101    double value;
     102
     103    snprintf (field, 256, "TS_%03d", i);
     104    if (!gfits_scan (phu, field, "%lf", 1, &value)) return FALSE;
     105    cam->tstart[i] = ohana_mjd_to_sec(value);
     106
     107    snprintf (field, 256, "TE_%03d", i);
     108    if (!gfits_scan (phu, field, "%lf", 1, &value)) return FALSE;
     109    cam->tstop[i] = ohana_mjd_to_sec(value);
     110  }
     111
    88112  while (TRUE) {
    89     Header *header = NULL;
    90     ALLOCATE (header, H
    91 
    92113    // load data for this header : if not found, assume we hit the end of the file
    93     if (!gfits_load_header (f, &theader)) break;
     114    Header *header = gfits_alloc_header();
     115    if (!gfits_load_header (f, header)) {
     116      free (header);
     117      break;
     118    }
    94119   
    95     if (!gfits_scan (&theader, "FILTER", "%d", 1, &filter)) return FALSE;
    96     if (!gfits_scan (&theader, "SEASON", "%d", 1, &season)) return FALSE;
    97 
     120    if (!gfits_scan (header, "FILTER", "%d", 1, &filter)) myAbort ("failed to find FILTER");
     121    if (!gfits_scan (header, "SEASON", "%d", 1, &season)) myAbort ("failed to find SEASON");
     122    if (!gfits_scan (header, "X_CHIP", "%d", 1, &ix))     myAbort ("failed to find X_CHIP");
     123    if (!gfits_scan (header, "Y_CHIP", "%d", 1, &iy))     myAbort ("failed to find Y_CHIP");
    98124    // XXX NOTE: astroflat.20150209.fits had ix and iy backwards in header
    99125    // double-check that the new flats are OK
    100    
    101     if (!gfits_scan (&theader, "X_CHIP",  "%d", 1, &ix))      return FALSE;
    102     if (!gfits_scan (&theader, "Y_CHIP",  "%d", 1, &iy))      return FALSE;
    103 
    104     Matrix *matrix = NULL;
    105     ALLOCATE (matrix, Matrix, 1);
    106     if (!gfits_load_matrix (f, matrix, &theader)) break;
     126
     127    Matrix *matrix = gfits_alloc_matrix();
     128    if (!gfits_load_matrix (f, matrix, header)) myAbort ("failed to read matrix");
    107129
    108130    int index = ix + iy*cam->Nx + filter*cam->Nchips + season*cam->Nflats;
     
    112134    myAssert (!cam->matrix[index], "entry already assigned?");
    113135
    114     // assert that cam->matrix[index] is NULL?
    115136    cam->matrix[index] = matrix;
    116137    cam->header[index] = header;
    117138  }
     139  return cam;
     140}
     141
     142int save_cam_correction (CamCorrection *cam, char *filename) {
     143
     144  int i;
     145
     146  /* open file for input */
     147  FILE *f = fopen (filename, "w");
     148  if (!f) {
     149    gprint (GP_ERR, "can't open file for write : %s\n", filename);
     150    return FALSE;
     151  }
     152
     153  Matrix matrix;
     154
     155  gfits_modify_alt (cam->phu, "EXTEND", "%t",  1, TRUE);
     156  gfits_fwrite_header  (f, cam->phu);
     157  gfits_create_matrix (cam->phu, &matrix);
     158  gfits_fwrite_matrix  (f, &matrix);
     159  // write an empty matrix?
     160
     161  for (i = 0; i < cam->Nvalues; i++) {
     162    if (!cam->matrix[i]) continue;
     163    gfits_primary_to_extended (cam->header[i], "IMAGE", "Image Extension");
     164    gfits_fwrite_header  (f, cam->header[i]);
     165    gfits_fwrite_matrix  (f, cam->matrix[i]);
     166    fflush (f);
     167  }
     168
    118169  return TRUE;
    119170}
    120171
    121 int save_cam_correction (char *filename) {
    122 
    123   int i, ix, iy, dir, filter;
    124   Header header;
    125 
    126   /* open file for input */
    127   FILE *f = fopen (filename, "r");
    128   if (f == (FILE *) NULL) {
    129     gprint (GP_ERR, "can't open file for read : %s\n", filename);
    130     return FALSE;
    131   }
    132 
    133   gfits_init_header (&header);
    134   header.extend = TRUE;
    135   gfits_create_header (&header);
    136   gfits_create_matrix (&header, &matrix);
    137 
    138   gfits_modify (&header, "NSEASON",  "%d", 1, camcoor->Nseason);
    139   gfits_modify (&header, "NFILTER",  "%d", 1, camcoor->Nfilter);
    140   gfits_modify (&header, "NX",       "%d", 1, camcoor->Nx);
    141   gfits_modify (&header, "NY",       "%d", 1, camcoor->Ny);
    142   gfits_modify (&header, "DX",       "%d", 1, camcoor->dX);
    143   gfits_modify (&header, "DY",       "%d", 1, camcoor->dY);
    144   gfits_modify (&header, "NX_CHIP",  "%d", 1, camcoor->NxCCD);
    145   gfits_modify (&header, "NY_CHIP",  "%d", 1, camcoor->NyCCD);
    146 
    147   gfits_fwrite_header  (f, &header);
    148   gfits_fwrite_matrix  (f, &matrix);
    149   gfits_free_header (&header);
    150   gfits_free_matrix (&matrix);
    151 
    152   for (i = 0; i < camcorr->Nvalues; i++) {
    153 
    154     // I don't need to do this, do I?
    155     // if (!gfits_scan (&theader, "FILTER", "%d", 1, &filter));
    156     // if (!gfits_scan (&theader, "SEASON", "%d", 1, &season));
    157     // if (!gfits_scan (&theader, "X_CHIP",  "%d", 1, &ix));   
    158     // if (!gfits_scan (&theader, "Y_CHIP",  "%d", 1, &iy));   
    159 
    160     if (!cam->matrix[i]) continue;
    161 
    162     gfits_fwrite_header  (f, whatHeader);
    163     gfits_fwrite_matrix  (f, cam->matrix[i]);
    164   }
    165 
    166   return TRUE;
    167 }
    168 
     172# if (0)
    169173// input is:
    170174// chipID == octal chip ID (eg 40 for XY40)
    171175// filter == (01234 = grizy)
    172 int get_cam_correction (int chipID, int filter, float Xccd, float Yccd, double *dX, double *dY) {
     176int get_cam_correction (CamCorrection *cam, int chipID, int filter, float Xccd, float Yccd, double *dX, double *dY) {
    173177
    174178  *dX = 0.0;
     
    215219  return TRUE;
    216220}
    217 
    218 float get_cam_correction_by_flatid (int flat_id, float Xccd, float Yccd) {
    219 
    220   float dM = 0.0;
     221# endif
     222
     223float get_cam_correction (CamCorrection *cam, int flat_id, float Xccd, float Yccd) {
     224
     225  if (!flat_id) return 0.0;
     226
     227  myAssert (flat_id > 0, "flat_id out of range");
     228  myAssert (flat_id <= cam->Nvalues, "flat_id out of range");
    221229
    222230  // validate the flat_id (not out of range?)
    223   int seq = camcorr->IDtoSeq[flat_id];
     231  // NOTE: we are setting the ID to be the sequence + 1 (match_camcorr_to_images.c:39)
     232  int seq = flat_id - 1;
    224233
    225234  Matrix *matrix = cam->matrix[seq];
    226   NxModel = cam->matrix[index]->Naxis[0];
    227 
    228   buffer = (float *) matrix->buffer;
    229 
    230   jx = MAX(MIN(cam->NxCCD, Xccd / cam->dX), 0);
    231   jy = MAX(MIN(cam->NyCCD, Yccd / cam->dY), 0);
    232 
    233   value = buffer[jx + jy*NxModel];
    234 
    235   dM = isnan(value) ? 0.0 : value;
     235  myAssert (matrix, "oops");
     236
     237  int NxModel = matrix->Naxis[0];
     238
     239  float *buffer = (float *) matrix->buffer;
     240
     241  int jx = MAX(MIN(cam->NxCCD, Xccd / cam->dX), 0);
     242  int jy = MAX(MIN(cam->NyCCD, Yccd / cam->dY), 0);
     243
     244  float value = buffer[jx + jy*NxModel];
     245
     246  float dM = isnan(value) ? 0.0 : value;
    236247
    237248  return dM;
    238249}
    239250
    240 int merge_flatcorr_and_camcorr (FlatCorrectionTable *flatcorr, CamCorrection *camcorr) {
     251CamCorrection *merge_flatcorr_and_camcorr (FlatCorrectionTable *flatcorr, CamCorrection *camcorr) {
     252
     253  int i, j, ix, iy, Nf, Ns;
    241254
    242255  // this function is specifically tuned for merging the ubercal flats with the high-res flats
     
    245258  // we have a highres flat with Nfilter, Nx, Ny images
    246259
    247   myAssert (flatcorr->Nfilter <= camcorr->Nfilter, "unclear how to use more flatcorr filters");
    248260  myAssert (camcorr->Nseason == 1, "no code to merge seasons flats to seasonal flats");
    249261
    250262  // we need to create a new CamCorrection with Nx,Ny,Nfilter from camcorr and Nseason from flatcorr:
    251263  CamCorrection *newcorr = alloc_cam_correction (camcorr->Nx, camcorr->Ny, camcorr->Nfilter, flatcorr->Nseason);
     264  newcorr->phu = gfits_alloc_header();
     265  gfits_create_header(newcorr->phu);
     266
    252267  newcorr->dX      = camcorr->dX;    // superpixel sampling, x-dir
    253268  newcorr->dY      = camcorr->dY;    // superpixel sampling, y-dir
     
    258273    newcorr->tstart[Ns] = flatcorr->tstart[Ns];
    259274    newcorr->tstop[Ns]  = flatcorr->tstop[Ns];
     275  }
     276
     277  gfits_modify (newcorr->phu, "NSEASON",  "%d", 1, newcorr->Nseason);
     278  gfits_modify (newcorr->phu, "NFILTER",  "%d", 1, newcorr->Nfilter);
     279  gfits_modify (newcorr->phu, "NX",       "%d", 1, newcorr->Nx);
     280  gfits_modify (newcorr->phu, "NY",       "%d", 1, newcorr->Ny);
     281  gfits_modify (newcorr->phu, "DX",       "%d", 1, newcorr->dX);
     282  gfits_modify (newcorr->phu, "DY",       "%d", 1, newcorr->dY);
     283  gfits_modify (newcorr->phu, "NX_CHIP",  "%d", 1, newcorr->NxCCD);
     284  gfits_modify (newcorr->phu, "NY_CHIP",  "%d", 1, newcorr->NyCCD);
     285
     286  // XXX UNITS for tstart, tstop
     287  for (Ns = 0; Ns < newcorr->Nseason; Ns++) {
     288    char field[256];
     289    double value;
     290    snprintf (field, 256, "TS_%03d", Ns);
     291    value = ohana_sec_to_mjd (newcorr->tstart[Ns]);
     292    if (!gfits_modify (newcorr->phu, field, "%lf", 1, value)) return FALSE;
     293
     294    snprintf (field, 256, "TE_%03d", Ns);
     295    value = ohana_sec_to_mjd (newcorr->tstop[Ns]);
     296    if (!gfits_modify (newcorr->phu, field, "%lf", 1, value))  return FALSE;
    260297  }
    261298
     
    265302      for (ix = 0; ix < newcorr->Nx; ix++) {
    266303        for (iy = 0; iy < newcorr->Ny; iy++) {
    267           int index = ix + iy*newcorr->Nx + Nf*newcorr->Nchips + Ns*newcorr->Nflats;
    268           ALLOCATE (newcorr->matrix[index], Matrix, 1);
    269          
    270304          // camcorr only had one season
     305          int idxNew = ix + iy*newcorr->Nx + Nf*newcorr->Nchips + Ns*newcorr->Nflats;
    271306          int idxOld = ix + iy*camcorr->Nx + Nf*camcorr->Nchips;
    272307
    273           // allocates the buffer
    274           gfits_copy_matrix (camcorr->matrix[idxOld], newcorr->matrix[index]);
     308          if (!camcorr->matrix[idxOld]) continue;
     309          myAssert (camcorr->header[idxOld], "oops");
     310
     311          // these allocate the structures but not the data
     312          newcorr->matrix[idxNew] = gfits_alloc_matrix();
     313          newcorr->header[idxNew] = gfits_alloc_header();
     314
     315          // these allocate the buffers
     316          gfits_copy_header (camcorr->header[idxOld], newcorr->header[idxNew]);
     317          gfits_copy_matrix (camcorr->matrix[idxOld], newcorr->matrix[idxNew]);
     318
     319          // using hard-wired filter names:
     320          if (!gfits_modify (newcorr->header[idxNew], "FILTNAME", "%s", 1, filters[Nf])) myAbort ("failed to set FILTNAME");
     321          if (!gfits_modify (newcorr->header[idxNew], "FILTER", "%d", 1, Nf)) myAbort ("failed to set FILTER");
     322          if (!gfits_modify (newcorr->header[idxNew], "SEASON", "%d", 1, Ns)) myAbort ("failed to set SEASON");
     323          if (!gfits_modify (newcorr->header[idxNew], "X_CHIP", "%d", 1, ix)) myAbort ("failed to set X_CHIP");
     324          if (!gfits_modify (newcorr->header[idxNew], "Y_CHIP", "%d", 1, iy)) myAbort ("failed to set Y_CHIP");
     325
     326          char extname[80];
     327          snprintf (extname, 80, "td_dM_%d_%s_%d_%d", Ns, filters[Ns], ix, iy);
     328          if (!gfits_modify (newcorr->header[idxNew], "EXTNAME", "%s", 1, extname)) myAbort ("failed to set EXTNAME");
    275329        }
    276330      }
     
    283337
    284338  // loop over the flatcorr images and match from image[seq] to newcorr:
    285 
    286339  for (i = 0; i < flatcorr->Nimage; i++) {
    287340
    288341    int ID = flatcorr->image[i].ID;
    289     myAssert (flatcorrTable->IDtoSeq[ID] == i, "oops");
     342    myAssert (flatcorr->IDtoSeq[ID] == i, "oops");
    290343
    291344    // find the season (these are inherited from flatcorr, so they must match
     
    300353
    301354    // find the chip and filter from photcode:
    302     iy = flatcorr->image[i].photcode % 10;
    303     ix = (int)(flatcorr->image[i].photcode / 10) % 10;
    304     filter = (int)(flatcorr->image[i].photcode / 100) % 10;
     355    int iy = photcode % 10;
     356    int ix = (int)(photcode / 10) % 10;
     357    int filter = (int)(photcode / 100) % 10;
    305358             
     359    myAssert (ix     < newcorr->Nx, "oops");
     360    myAssert (iy     < newcorr->Ny, "oops");
     361    myAssert (filter < newcorr->Nfilter, "oops");
     362
    306363    // seq is the flatcorr entry of interest
    307     int seq = ix + iy * newcorr->Nx + filter * newcorr->Nchips + season * camcorr->Nflat;
     364    int seq = ix + iy * newcorr->Nx + filter * newcorr->Nchips + season * newcorr->Nflats;
    308365
    309366    // we now have newcorr->matrix[seq]
     367    // flatcorr contains entries for unpopulated corner chips, but camcorr does not
    310368    Matrix *newDelta = newcorr->matrix[seq];
     369    if (!newDelta) continue;
     370
     371    float *values = (float *) newDelta->buffer;
    311372
    312373    int jx, jy;
     374    int NxModel = newDelta->Naxis[0];
     375    int NyModel = newDelta->Naxis[1];
    313376
    314377    // now loop over the newcorr pixels to get the flatcorr value:
    315     for (jx = 0; jx < newcorr->NxCCD; jx++) {
    316       for (jy = 0; jy < newcorr->NyCCD; jy++) {
     378    for (jx = 0; jx < NxModel; jx++) {
     379      for (jy = 0; jy < NyModel; jy++) {
    317380       
    318381        int chipX = jx * newcorr->dX;
     
    320383
    321384        // get the flatcorr offset from this chip coordinate
    322         float delta = FlatCorrectionOffset (flatcorr, flatcorr->image[i].ID, chipX, chipY);
    323 
    324         newDelta->buffer[jx + jy*NxCCD] += flatcorr->offset[i][Jx][Jy];
     385        float delta = FlatCorrectionOffset (flatcorr, ID, chipX, chipY);
     386
     387        values[jx + jy*NxModel] += delta;
    325388      }
    326389    }
    327390  }   
    328391
    329   return TRUE;
    330 }
    331 
    332 
     392  return newcorr;
     393}
     394
     395
  • trunk/Ohana/src/uniphot/src/initialize_setphot.c

    r39262 r39272  
    3535
    3636  CAM_FILE = NULL;
    37   if ((N = get_argument (argc, argv, "-CAM"))) {
     37  if ((N = get_argument (argc, argv, "-cam-flat"))) {
    3838    remove_argument (N, &argc, argv);
    3939    char *tmpfile = strcreate (argv[N]);
    4040    CAM_FILE = abspath (tmpfile, DVO_MAX_PATH);
    4141    remove_argument (N, &argc, argv);
    42   }
    43   CAM_RESET = FALSE;
    44   if ((N = get_argument (argc, argv, "-CAM-reset"))) {
    45     remove_argument (N, &argc, argv);
    46     CAM_RESET = TRUE;
    4742  }
    4843
     
    131126    remove_argument (N, &argc, argv);
    132127    UBERCAL = TRUE;
     128    if (!CAM_FILE) {
     129      fprintf (stderr, "ubercal now requires a high-res static flat-file: -cam-flat filename\n");
     130      exit (2);
     131    }
    133132  }
    134133
  • trunk/Ohana/src/uniphot/src/load_zpt_table.c

    r39138 r39272  
    253253    } else {
    254254      if (!gfits_scan (&header, "FILTER", "%s", 1, filter)) {
    255         fprintf (stderr, "cannot find FILTER in header of %s, using %s\n", filename, filters_uc[nfilter]);
     255        fprintf (stderr, "warning: cannot find FILTER in header of %s, using %s\n", filename, filters_uc[nfilter]);
    256256        strcpy (filter, filters_uc[nfilter]);
    257257        // fclose (f);
  • trunk/Ohana/src/uniphot/src/match_camcorr_to_images.c

    r39263 r39272  
    44// the photcode is used to find the actual flat correction
    55
    6 int match_camcorr_to_images (Image *image, off_t Nimage) {
     6int match_camcorr_to_images (Image *image, off_t Nimage, CamCorrection *camcorr) {
    77
    88  int i, j;
    9 
    10   CamCorrection *camcorr = get_cam_correction_ptr (); 
    119
    1210  // we have an array of CamCorrection->matrix values, where each matrix is a flat-field image
     
    1917  int maxCode = 10576;
    2018
     19  int Nmissed = 0;
     20 
    2121  for (i = 0; i < Nimage; i++) {
    2222    if (!image[i].photcode) continue; // skip PHU images
     
    3030      if (image[i].tzero > camcorr->tstop[j]) continue;
    3131
    32       int iy = image[i].photcode % 10;
    33       int ix = (int)(image[i].photcode / 10) % 10;
    34       int filter = (int)(image[i].photcode / 100) % 10;
     32      int photcode = image[i].photcode;
     33      int iy =  photcode % 10;
     34      int ix = (int)(photcode / 10) % 10;
     35      int filter = (int)(photcode / 100) % 10;
    3536
    36       int seq = ix + iy * camcorr->Nx + filter * camcorr->Nchips + j * camcorr->Ngroup;
    37 
    38       image[i].photom_map_id = camcorr->imageID[seq].ID;
     37      // we add one so photom_map_id = 0 can mean no map
     38      int seq = ix + iy * camcorr->Nx + filter * camcorr->Nchips + j * camcorr->Nflats;
     39      image[i].photom_map_id = seq + 1;
    3940      found = TRUE;
    4041    }
  • trunk/Ohana/src/uniphot/src/setphot.c

    r39267 r39272  
    77  FITS_DB db;
    88  ZptTable *zpts;
    9   Image *image;
    109
    11   // FlatCorrectionTable flatcorrTable;
     10  CamCorrection *camcorr = NULL;
     11  FlatCorrectionTable flatcorrTable;
    1212
    1313  /* get configuration info, args, lockfile */
     
    3131    // XXX should this program save any old copy of this file?
    3232
    33     if (!load_cam_correction (CAM_FILE)) {
     33    CamCorrection *rawcorr = load_cam_correction (CAM_FILE);
     34    if (!rawcorr) {
    3435      fprintf (stderr, "failed to load camera-static flat-field correction\n");
    3536      exit (2);
    3637    }
    37     merge_flatcorr_with_cam (flatcorrTable);
     38    camcorr = merge_flatcorr_and_camcorr (&flatcorrTable, rawcorr);
    3839
    3940    char newflatfile[DVO_MAX_PATH];
    40     int size = snprintf (newflatfile, DVO_MAX_PATH, "%s/flatfield.fits", CATDIR);
     41    size = snprintf (newflatfile, DVO_MAX_PATH, "%s/flatfield.fits", CATDIR);
    4142    assert (size < DVO_MAX_PATH);
    42     save_cam_correction (newflatfile);
     43    save_cam_correction (camcorr, newflatfile);
    4344  } else {
    4445    zpts = load_zpt_table (argv[1], &Nzpts);
     
    4849
    4950  // load images
    50   image  = load_images_setphot (&db, &Nimage);
     51  Image *image  = load_images_setphot (&db, &Nimage);
    5152  if (!UPDATE) dvo_image_unlock (&db);
    5253 
     
    5657    // we are going to deprecate the flatcorr imaage lookup
    5758    // match_flatcorr_to_images (image, Nimage, &flatcorrTable);
    58     match_camcorr_to_images (image, Nimage);
     59    match_camcorr_to_images (image, Nimage, camcorr);
    5960  }
    6061
    6162  if (!IMAGES_ONLY) {
    62     status = update_dvo_setphot (image, Nimage);
     63    status = update_dvo_setphot (image, Nimage, camcorr);
    6364  }
    6465
  • trunk/Ohana/src/uniphot/src/update_catalog_setphot.c

    r39262 r39272  
    11# include "setphot.h"
    22
    3 void update_catalog_setphot (Catalog *catalog, Image *image, off_t *index, off_t Nimage) {
     3void update_catalog_setphot (Catalog *catalog, Image *image, off_t *index, off_t Nimage, CamCorrection *camcorr) {
    44
    55  off_t i, j, found;
     
    77  // if we are resetting, reset all flags
    88  DVOMeasureFlags photomFlags = ID_MEAS_POOR_PHOTOM | ID_MEAS_SKIP_PHOTOM | ID_MEAS_PHOTOM_UBERCAL;
    9 
    10   CamCorrection *camcorr = get_cam_correction_ptr (); 
    119
    1210  found = 0;   
     
    2826      int flat_id = image[id].photom_map_id;
    2927      if (flat_id > 0) {
    30         Mflat = CamCorrectionOffset (camcorr, flat_id, catalog[0].measure[m].Xccd, catalog[0].measure[m].Yccd);
     28        Mflat = get_cam_correction (camcorr, flat_id, catalog[0].measure[m].Xccd, catalog[0].measure[m].Yccd);
    3129      }
    3230
  • trunk/Ohana/src/uniphot/src/update_dvo_setphot.c

    r38986 r39272  
    33
    44// XXX convert this to use the ImageSubset to save space
    5 int update_dvo_setphot (Image *image, off_t Nimage, FlatCorrectionTable *flatcorr) {
     5int update_dvo_setphot (Image *image, off_t Nimage, CamCorrection *camcorr) {
    66
    77  SkyTable *sky = NULL;
     
    7373    }
    7474
    75     update_catalog_setphot (&catalog, image, index, Nimage, flatcorr);
     75    update_catalog_setphot (&catalog, image, index, Nimage, camcorr);
    7676
    7777    if (!UPDATE) {
Note: See TracChangeset for help on using the changeset viewer.