IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 37913


Ignore:
Timestamp:
Feb 16, 2015, 8:17:35 PM (11 years ago)
Author:
eugene
Message:

adding cam correction options

Location:
branches/eam_branches/ipp-20150112/Ohana/src/uniphot
Files:
6 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20150112/Ohana/src/uniphot/Makefile

    r37807 r37913  
    112112$(SRC)/setastrom.$(ARCH).o          \
    113113$(SRC)/initialize_setastrom.$(ARCH).o \
     114$(SRC)/cam_correction.$(ARCH).o \
    114115$(SRC)/dcr_correction.$(ARCH).o \
    115116$(SRC)/kh_correction.$(ARCH).o \
    … …  
    127128$(SRC)/update_dvo_setastrom.$(ARCH).o \
    128129$(SRC)/update_catalog_setastrom.$(ARCH).o \
     130$(SRC)/cam_correction.$(ARCH).o \
    129131$(SRC)/dcr_correction.$(ARCH).o \
    130132$(SRC)/kh_correction.$(ARCH).o \

  • branches/eam_branches/ipp-20150112/Ohana/src/uniphot/include/setastrom.h

    r37910 r37913  
    1313// we have one correction (an image) for each filter and chip
    1414typedef struct {
    15   int Nvalues;
    16   int Nfilters; // number of filters
    17   int NchipX;   // number of chips in x
    18   int NchipY;   // number of chips in y
     15  int Nx;       // number of chips in x
     16  int Ny;       // number of chips in y
     17  int Nfilter;  // number of filters
     18  int Ndir;     // number of correction dimensions
     19
     20  int Nchips;   // chip offset (Nx*Ny)
     21  int Ngroup;   // direction offset (Nx*Ny*Nfilters)
     22  int Nvalues;  // Nx*Ny*Nfilters*Ndir
     23
     24  int dX;       // superpixel size
     25  int dY;       // superpixel size
     26
    1927  int NxCCD;    // number of pixels
     28  int NyCCD;    // number of pixels
     29
    2030  Matrix **matrix; // allocate an array of pointers
    21   // index = ix + iy*Nx + filter*Nx*Ny + dir*Nx*Ny*Nfilter
     31  // index = ix + iy*Nx + filter*Nchips + dir*Ngroup
    2232} CamCorrection;
    2333
    … …  
    3848char        *KH_FILE;
    3949char        *DCR_FILE;
     50char        *CAM_FILE;
    4051
    4152int          KH_RESET;
    4253int          DCR_RESET;
     54int          CAM_RESET;
    4355
    4456SkyRegion    UserPatch;
    … …  
    7688int load_dcr_correction (char *filename);
    7789int get_dcr_correction (int filter, double *dX, double *dY, float Color);
     90
     91int load_cam_correction (char *filename);
     92int get_cam_correction (int chipID, int filter, float Xccd, float Yccd, double *dX, double *dY);

  • branches/eam_branches/ipp-20150112/Ohana/src/uniphot/src/cam_correction.c

    r37910 r37913  
    66// note that we are allocating pointers for XY00 - XY77, but only the octal elements are set (and not the 00,07,70,77 ones)
    77
    8 CamCorrection *load_cam_correction (char *filename) {
     8static CamCorrection *cam = NULL;
    99
    10   int Ncol;
    11   off_t Nrow;
     10int load_cam_correction (char *filename) {
    1211
    13   char type[16];
    14   char extname[80];
    15 
     12  int i, ix, iy, dir, filter;
    1613  Header header;
    17   Header theader;
    18   FTable ftable;
    1914
    2015  /* open file for input */
    … …  
    2924
    3025  // read the PHU header
    31   if (!gfits_load_header (f, &header)) return (FALSE);
     26  if (!gfits_load_header (f, &header)) return FALSE;
    3227 
    3328  int Nbytes = gfits_data_size (&header);
    3429  fseeko (f, Nbytes, SEEK_CUR);
    3530
    36   ftable.header = &theader;
     31  int Ndir, Nfilter, Nx, Ny, dX, dY, NxCCD, NyCCD;
     32  if (!gfits_scan (&header, "NDIR",    "%d", 1, &Ndir))    return FALSE;
     33  if (!gfits_scan (&header, "NFILTER", "%d", 1, &Nfilter)) return FALSE;
     34  if (!gfits_scan (&header, "NX",       "%d", 1, &Nx))      return FALSE;
     35  if (!gfits_scan (&header, "NY",       "%d", 1, &Ny))      return FALSE;
     36  if (!gfits_scan (&header, "DX",       "%d", 1, &dX))      return FALSE;
     37  if (!gfits_scan (&header, "DY",       "%d", 1, &dY))      return FALSE;
    3738
    38   if (!gfits_scan (&theader, "NDIR",    "%d", 1, &Ndir))    return (NULL);
    39   if (!gfits_scan (&theader, "NFILTER", "%d", 1, &Nfilter)) return (NULL);
    40   if (!gfits_scan (&theader, "NX_CHIP", "%d", 1, &Nx))      return (NULL);
    41   if (!gfits_scan (&theader, "NY_CHIP", "%d", 1, &Ny))      return (NULL);
    42   if (!gfits_scan (&theader, "DX_CHIP", "%d", 1, &dX))      return (NULL);
    43   if (!gfits_scan (&theader, "DY_CHIP", "%d", 1, &dY))      return (NULL);
     39  if (!gfits_scan (&header, "NX_CHIP", "%d", 1, &NxCCD))   return FALSE;
     40  if (!gfits_scan (&header, "NY_CHIP", "%d", 1, &NyCCD))   return FALSE;
    4441
    4542  ALLOCATE (cam->matrix, Matrix *, Ndir*Nfilter*Nx*Ny);
    … …  
    4744  cam->Nx      = Nx;       // for gpc1, should be 8
    4845  cam->Ny      = Ny;       // for gpc1, should be 8
     46  cam->Nfilter = Nfilter;  // for gpc1, should be 5
     47  cam->Ndir    = Ndir;     // for gpc1, should be 2
     48
    4949  cam->dX      = dX;       // superpixel sampling, x-dir
    5050  cam->dY      = dY;       // superpixel sampling, y-dir
    51   cam->Ndir    = Ndir;     // for gpc1, should be 2
    52   cam->Nfilter = Nfilter;  // for gpc1, should be 5
    5351
     52  cam->NxCCD   = NxCCD;    // pixels per chip
     53  cam->NyCCD   = NyCCD;    // pixels per chip
     54
     55  cam->Nchips  = Nx * Ny;
     56  cam->Ngroup  = Nx * Ny * Nfilter;
     57  cam->Nvalues = Nx * Ny * Nfilter * Ndir;
     58
     59  ALLOCATE (cam->matrix, Matrix *, cam->Nvalues);
     60  for (i = 0; i < cam->Nvalues; i++) {
     61    cam->matrix[i] = NULL;
     62  }
     63
     64  int found = FALSE;
    5465  while (TRUE) {
     66    Header theader;
     67
    5568    // load data for this header : if not found, assume we hit the end of the file
    5669    if (!gfits_load_header (f, &theader)) break;
    5770   
    58     if (!gfits_scan (&theader, "FILTER",  "%d", 1, &filter))  return (FALSE);
    59     if (!gfits_scan (&theader, "DIR",     "%d", 1, &dir))     return (FALSE);
    60     if (!gfits_scan (&theader, "X_CHIP",  "%d", 1, &ix))      return (FALSE);
    61     if (!gfits_scan (&theader, "Y_CHIP",  "%d", 1, &iy))      return (FALSE);
     71    if (!gfits_scan (&theader, "FILTER",  "%d", 1, &filter))  return FALSE;
     72    if (!gfits_scan (&theader, "DIR",     "%d", 1, &dir))     return FALSE;
     73    if (!gfits_scan (&theader, "X_CHIP",  "%d", 1, &ix))      return FALSE;
     74    if (!gfits_scan (&theader, "Y_CHIP",  "%d", 1, &iy))      return FALSE;
    6275
    6376    Matrix *matrix = NULL;
    6477    ALLOCATE (matrix, Matrix, 1);
    65     if (!gfits_load_header (f, matrix, &theader)) break;
     78    if (!gfits_load_matrix (f, matrix, &theader)) break;
    6679
    67     int index = ix + iy*Nx + filter*Nx*Ny + dir*Nx*Ny*Nfilter;
     80    int index = ix + iy*cam->Nx + filter*cam->Nchips + dir*cam->Ngroup;
    6881    myAssert (index >= 0, "index too small");
    69     myAssert (index < Ndir*Nfilter*Nx*Ny, "index too big");
     82    myAssert (index < cam->Nvalues, "index too big");
     83
     84    myAssert (!cam->matrix[i], "entry already assigned?");
    7085
    7186    // assert that cam->matrix[index] is NULL?
    7287    cam->matrix[index] = matrix;
     88
     89    if (!found) {
     90     
     91    }
    7392  }
    74   return (TRUE);
     93  return TRUE;
    7594}
    7695
    77 int get_cam_correction (CamCorrection *cam, int chip, int filter, double *dX, double *dY, float Xccd, float Yccd) {
     96// input is:
     97// chipID == octal chip ID (eg 40 for XY40)
     98// filter == (01234 = grizy)
     99int get_cam_correction (int chipID, int filter, float Xccd, float Yccd, double *dX, double *dY) {
    78100
    79101  *dX = 0.0;
    80102  *dY = 0.0;
    81103
    82   // split out the chip number (eg 43 for XY43) into X and Y elements
    83   int ix = (int) (chip / 10);
    84   int iy = (int) (chip % 10);
     104  int index, jx, jy, NxModel;
     105  Matrix *matrix;
     106  float *buffer, value;
    85107
    86   int dir = 0; // dX
    87   int index = ix + iy*cam->Nx + filter*cam->Nx*cam->Ny + dir*cam->Nx*cam->Ny*cam->Nfilter;
     108  // split out the chipID number (eg 43 for XY43) into X and Y elements
     109  int ix = (int) (chipID / 10);
     110  int iy = (int) (chipID % 10);
    88111
    89   Matrix *matrix = cam->matrix[index];
     112  // dX (dir == 0)
     113  index = ix + iy*cam->Nx + filter*cam->Nchips;
    90114
    91   float *buffer = (float *) matrix->buffer;
     115  matrix = cam->matrix[index];
     116  NxModel = cam->matrix[index]->Naxis[0];
    92117
    93   int jx = MAX(MIN(Nx, Xccd / cam->dX), 0);
    94   int jy = MAX(MIN(Nx, Yccd / cam->dY), 0);
     118  buffer = (float *) matrix->buffer;
    95119
    96   float value = buffer[jx + jy*Nx];
     120  jx = MAX(MIN(cam->NxCCD, Xccd / cam->dX), 0);
     121  jy = MAX(MIN(cam->NyCCD, Yccd / cam->dY), 0);
     122
     123  value = buffer[jx + jy*NxModel];
    97124
    98125  *dX = isnan(value) ? 0.0 : value;
    99   return (TRUE);
     126
     127  // dY (dir == 1) uses the next block
     128  index += cam->Ngroup;
     129
     130  matrix = cam->matrix[index];
     131  NxModel = cam->matrix[index]->Naxis[0];
     132
     133  buffer = (float *) matrix->buffer;
     134
     135  jx = MAX(MIN(cam->NxCCD, Xccd / cam->dX), 0);
     136  jy = MAX(MIN(cam->NyCCD, Yccd / cam->dY), 0);
     137
     138  value = buffer[jx + jy*NxModel];
     139
     140  *dY = isnan(value) ? 0.0 : value;
     141
     142  return TRUE;
    100143}
    101 

  • branches/eam_branches/ipp-20150112/Ohana/src/uniphot/src/initialize_setastrom.c

    r37807 r37913  
    44    fprintf (stderr, "USAGE: setastrom [options]\n");
    55    fprintf (stderr, "  options:\n");
     6    fprintf (stderr, "    -KH (file) : supply Koppenhoefer correction splines\n");
    67    fprintf (stderr, "    -DCR (file) : supply DCR correction splines\n");
    7     fprintf (stderr, "    -KH (file) : supply Koppenhoefer correction splines\n");
     8    fprintf (stderr, "    -CAM (file) : supply camera-static correction file\n");
    89    fprintf (stderr, "    -v : verbose mode\n");
    910    fprintf (stderr, "    -region Rmin Rmax Dmin Dmax\n");
    … …  
    7778  }
    7879
     80  CAM_FILE = NULL;
     81  if ((N = get_argument (argc, argv, "-CAM"))) {
     82    remove_argument (N, &argc, argv);
     83    char *tmpfile = strcreate (argv[N]);
     84    CAM_FILE = abspath (tmpfile, DVO_MAX_PATH);
     85    remove_argument (N, &argc, argv);
     86  }
     87  CAM_RESET = FALSE;
     88  if ((N = get_argument (argc, argv, "-CAM-reset"))) {
     89    remove_argument (N, &argc, argv);
     90    CAM_RESET = TRUE;
     91  }
     92
    7993  SINGLE_CPT = NULL;
    8094  if ((N = get_argument (argc, argv, "-cpt"))) {
    … …  
    146160  }
    147161
    148   if (!KH_FILE && !DCR_FILE) {
    149     fprintf (stderr, "at least one of -KH and -DCR must be supplied\n");
    150     exit (1);
     162  if (!KH_FILE && !DCR_FILE && !CAM_FILE) {
     163    fprintf (stderr, "WARNING: none of -CAM, -KH, -DCR supplied\n");
    151164  }
    152165
    … …  
    159172  fprintf (stderr, "USAGE: setastrom_client -hostID (hostID) -catdir (catdir) -hostdir (hostdir) [options]\n");
    160173  fprintf (stderr, "  options:\n");
     174  fprintf (stderr, "    -KH  (file) : supply Koppenhoefer correction splines\n");
    161175  fprintf (stderr, "    -DCR (file) : supply DCR correction splines\n");
    162   fprintf (stderr, "    -KH (file) : supply Koppenhoefer correction splines\n");
     176  fprintf (stderr, "    -CAM (file) : supply camera-static correction file\n");
    163177  fprintf (stderr, "    -region Rmin Rmax Dmin Dmax\n");
    164178  fprintf (stderr, "    -update-catformat (format) : change database schema on output\n");
    … …  
    221235    remove_argument (N, &argc, argv);
    222236    DCR_RESET = TRUE;
     237  }
     238
     239  CAM_FILE = NULL;
     240  if ((N = get_argument (argc, argv, "-CAM"))) {
     241    remove_argument (N, &argc, argv);
     242    CAM_FILE = strcreate (argv[N]);
     243    remove_argument (N, &argc, argv);
     244  }
     245  CAM_RESET = FALSE;
     246  if ((N = get_argument (argc, argv, "-CAM-reset"))) {
     247    remove_argument (N, &argc, argv);
     248    CAM_RESET = TRUE;
    223249  }
    224250
    … …  
    282308  if (!CATDIR) usage_setastrom_client();
    283309
    284   if (!KH_FILE && !DCR_FILE) {
    285     fprintf (stderr, "at least one of -KH and -DCR must be supplied\n");
    286     exit (1);
     310  if (!KH_FILE && !DCR_FILE && !CAM_FILE) {
     311    fprintf (stderr, "WARNING: none of -CAM, -KH, -DCR supplied\n");
    287312  }
    288313

  • branches/eam_branches/ipp-20150112/Ohana/src/uniphot/src/update_catalog_setastrom.c

    r37908 r37913  
    121121    double dX_CAM = 0.0;
    122122    double dY_CAM = 0.0;
    123     if (CAM_FILE && doCAM) {
     123    if (CAM_FILE) {
    124124
    125125      // camera systematic correction ("astroflat") depends on only the X,Y coordinate,
    126126      // the filter, and the chip
    127127     
    128       int filtCode = (int)(measure->photcode / 100); // eg, 101 = r
     128      // 10134 : r, XY34 -> filtCode = 1
     129      int filtCode = (int)((measure->photcode % 100) / 100);
    129130
    130131      // correction is in arcseconds
    131       get_cam_correction (0, chipID, filtCode, measure->Xccd, measure->Yccd, &dX_CAM, &dY_CAM);
     132      get_cam_correction (chipID, filtCode, measure->Xccd, measure->Yccd, &dX_CAM, &dY_CAM);
    132133    }
    133134
    … …  
    260261    }
    261262    if (DCR_RESET) {
    262       measure[0].XoffDCR = -dX_DCR / pltScale;
    263       measure[0].YoffDCR = -dY_DCR / pltScale;
     263      measure[0].XoffDCR = 0.0;
     264      measure[0].YoffDCR = 0.0;
     265    }
     266    if (CAM_FILE) {
     267      measure[0].XoffCAM = -dX_DCR / pltScale;
     268      measure[0].YoffCAM = -dY_DCR / pltScale;
     269    }
     270    if (CAM_RESET) {
     271      measure[0].XoffDCR = 0.0;
     272      measure[0].YoffDCR = 0.0;
    264273    }
    265274
    … …  
    274283      measure[0].Yfix += measure[0].YoffDCR;
    275284    }
     285    if (isfinite(measure[0].XoffCAM) && isfinite(measure[0].YoffCAM)) {
     286      measure[0].Xfix += measure[0].XoffCAM;
     287      measure[0].Yfix += measure[0].YoffCAM;
     288    }
    276289    found ++;
    277290  }

  • branches/eam_branches/ipp-20150112/Ohana/src/uniphot/src/update_dvo_setastrom.c

    r37807 r37913  
    1717  }
    1818
    19   if (KH_FILE)  load_kh_correction (KH_FILE);
    20   if (DCR_FILE) load_dcr_correction (DCR_FILE);
     19  if (KH_FILE)  if (!load_kh_correction  (KH_FILE))  { fprintf (stderr, "failed to load KH  correction %s\n", KH_FILE ); exit (1); }
     20  if (DCR_FILE) if (!load_dcr_correction (DCR_FILE)) { fprintf (stderr, "failed to load DCR correction %s\n", DCR_FILE); exit (1); }
     21  if (CAM_FILE) if (!load_cam_correction (CAM_FILE)) { fprintf (stderr, "failed to load CAM correction %s\n", CAM_FILE); exit (1); }
    2122
    2223  // determine the populated SkyRegions overlapping the requested area (default depth)
    … …  
    104105    if (KH_FILE)          { snprintf (tmpline, 1024, "%s -KH %s",               command, KH_FILE);          strcpy (command, tmpline); }
    105106    if (DCR_FILE)         { snprintf (tmpline, 1024, "%s -DCR %s",              command, DCR_FILE);         strcpy (command, tmpline); }
     107    if (CAM_FILE)         { snprintf (tmpline, 1024, "%s -CAM %s",              command, CAM_FILE);         strcpy (command, tmpline); }
    106108    if (UPDATE_CATFORMAT) { snprintf (tmpline, 1024, "%s -update-catformat %s", command, UPDATE_CATFORMAT); strcpy (command, tmpline); }
    107109
    108110    if (KH_RESET)         { snprintf (tmpline, 1024, "%s -KH-reset",            command);                   strcpy (command, tmpline); }
    109111    if (DCR_RESET)        { snprintf (tmpline, 1024, "%s -DCR-reset",           command);                   strcpy (command, tmpline); }
     112    if (CAM_RESET)        { snprintf (tmpline, 1024, "%s -CAM-reset",           command);                   strcpy (command, tmpline); }
    110113
    111114    fprintf (stderr, "command: %s\n", command);
Note: See TracChangeset for help on using the changeset viewer.