IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 34844 for trunk


Ignore:
Timestamp:
Dec 18, 2012, 10:04:35 AM (14 years ago)
Author:
eugene
Message:

merging changes from eam_branches/ipp-20121130

Location:
trunk/Ohana
Files:
28 edited
2 copied

Legend:

Unmodified
Added
Removed

  • trunk/Ohana

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

    r34291 r34844  
    1515
    1616void usage (void) {
    17   fprintf (stderr, "USAGE: findcell -mktree (tree) (catdir)\n");
     17  fprintf (stderr, "USAGE: findcell -mktree (tree) (catdir) [-nx Nx] [-ny Ny]\n");
    1818  fprintf (stderr, "USAGE: findcell -tree (tree) (datafile)\n");
    1919  fprintf (stderr, "   (datafile) should contain a list of RA,DEC pairs\n");
    … …  
    2424int apply_tree (char *treefile, char *datafile);
    2525
     26float SCALE = 1.0;
     27int NX_SUB = 1;
     28int NY_SUB = 1;
     29
    2630int main (int argc, char **argv) {
    2731
    … …  
    3741  if (get_argument (argc, argv, "-help")) usage ();
    3842  if (get_argument (argc, argv, "-h")) usage ();
     43
     44  NX_SUB = NY_SUB = 1;
     45  if ((N = get_argument (argc, argv, "-nx"))) {
     46    remove_argument (N, &argc, argv);
     47    NX_SUB = atof (argv[N]);
     48    remove_argument (N, &argc, argv);
     49  }
     50  if ((N = get_argument (argc, argv, "-ny"))) {
     51    remove_argument (N, &argc, argv);
     52    NY_SUB = atof (argv[N]);
     53    remove_argument (N, &argc, argv);
     54  }
     55
     56  /* pixel scale (arcsec/pixel) */
     57  SCALE = 1.0;
     58  if ((N = get_argument (argc, argv, "-scale"))) {
     59    remove_argument (N, &argc, argv);
     60    SCALE = atof (argv[N]);
     61    remove_argument (N, &argc, argv);
     62  }
    3963
    4064  /* extra error messages */
    … …  
    109133
    110134  tree.Nzone = 46;
     135  tree.NX_SUB = NX_SUB;
     136  tree.NY_SUB = NY_SUB;
     137  tree.dPix = SCALE/3600.0;
    111138
    112139  ALLOCATE (tree.Nband, int, tree.Nzone);
    … …  
    122149  ALLOCATE (tree.ra,   double *, tree.Nzone);
    123150  ALLOCATE (tree.dec,  double *, tree.Nzone);
     151  ALLOCATE (tree.Xo,   double *, tree.Nzone);
     152  ALLOCATE (tree.Yo,   double *, tree.Nzone);
     153  ALLOCATE (tree.dX,      int *, tree.Nzone);
     154  ALLOCATE (tree.dY,      int *, tree.Nzone);
    124155  ALLOCATE (tree.cell,    int *, tree.Nzone);
    125156  ALLOCATE (tree.name,  char **, tree.Nzone);
    … …  
    152183    ALLOCATE (tree.ra[zone],   double, tree.NBAND[zone]);
    153184    ALLOCATE (tree.dec[zone],  double, tree.NBAND[zone]);
     185    ALLOCATE (tree.Xo[zone],   double, tree.NBAND[zone]);
     186    ALLOCATE (tree.Yo[zone],   double, tree.NBAND[zone]);
     187    ALLOCATE (tree.dX[zone],      int, tree.NBAND[zone]);
     188    ALLOCATE (tree.dY[zone],      int, tree.NBAND[zone]);
    154189    ALLOCATE (tree.cell[zone], int,    tree.NBAND[zone]);
    155190    ALLOCATE (tree.name[zone], char *, tree.NBAND[zone]);
    … …  
    157192      tree.ra[zone][band] = NAN;
    158193      tree.dec[zone][band] = NAN;
     194      tree.Xo[zone][band] = NAN;
     195      tree.Yo[zone][band] = NAN;
     196      tree.dX[zone][band] = -1;
     197      tree.dY[zone][band] = -1;
    159198      tree.cell[zone][band] = -1;
    160199      ALLOCATE (tree.name[zone][band], char, BOUNDARY_TREE_NAME_LENGTH);
    … …  
    179218      REALLOCATE (tree.ra[zone],   double, tree.NBAND[zone]);
    180219      REALLOCATE (tree.dec[zone],  double, tree.NBAND[zone]);
    181       REALLOCATE (tree.cell[zone], int,    tree.NBAND[zone]);
     220      REALLOCATE (tree.Xo[zone],   double, tree.NBAND[zone]);
     221      REALLOCATE (tree.Yo[zone],   double, tree.NBAND[zone]);
     222      REALLOCATE (tree.dX[zone],      int, tree.NBAND[zone]);
     223      REALLOCATE (tree.dY[zone],      int, tree.NBAND[zone]);
     224      REALLOCATE (tree.cell[zone],    int, tree.NBAND[zone]);
    182225      REALLOCATE (tree.name[zone], char *, tree.NBAND[zone]);
    183226      for (j = start; j < tree.NBAND[zone]; j++) {
    184227        tree.ra[zone][j] = NAN;
    185228        tree.dec[zone][j] = NAN;
     229        tree.Xo[zone][band] = NAN;
     230        tree.Yo[zone][band] = NAN;
     231        tree.dX[zone][band] = -1;
     232        tree.dY[zone][band] = -1;
    186233        tree.cell[zone][j] = -1;
    187234        ALLOCATE (tree.name[zone][j], char, BOUNDARY_TREE_NAME_LENGTH);
    … …  
    190237    tree.ra[zone][band] = ra;
    191238    tree.dec[zone][band] = dec;
     239    tree.Xo[zone][band] = x;
     240    tree.Yo[zone][band] = y;
     241    tree.dX[zone][band] = image[i].NX / NX_SUB;
     242    tree.dY[zone][band] = image[i].NY / NY_SUB;
     243   
    192244    tree.cell[zone][band] = i;
    193245
    … …  
    325377    }
    326378
    327     fprintf (stdout, "%10.6f %10.6f  %3d %3d  %s\n", ra, dec, zone, band, tree->name[zone][band]);
     379    // I know the projection cell (band,zone), but I need to find the skycell within that projection cell.
     380    // the proj cell is divided into Nx, Ny bits.
     381    // (ra,dec) for (Ro,Do) -> (x,y).  given (Xo,Yo),(dX,dY) I can find ix,iy
     382    // I currently track Ro,Do (tree->ra[zone][band], tree->dec[zone][band])
     383    // I need to have the skycell center in pixels (Xo,Yo) and the scale dX,dY
     384   
     385    // convert R,D to X,Y with hard-wired projection and scale, orientation?
     386   
     387    double x = 0.0;
     388    double y = 0.0;
     389    BoundaryTreeProjection (&x, &y, ra, dec, tree, zone, band);
     390
     391    int xi = x / tree->dX[zone][band];
     392    int yi = y / tree->dY[zone][band];
     393   
     394    char format[24], skycellname[128];
     395    int Ndigit = (int)(log10(tree->NX_SUB*tree->NY_SUB)) + 1 ;
     396    snprintf (format, 24, "%s.%%0%dd", tree->name[zone][band], Ndigit);
     397
     398    int N = xi + tree->NX_SUB * yi;
     399    snprintf (skycellname, 128, format, N);
     400
     401    fprintf (stdout, "%10.6f %10.6f  %8.3f %8.3f  %3d %3d  %s\n", ra, dec, x, y, zone, band, skycellname);
    328402  }
    329403
    330404  exit (0);
    331405}
     406

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

    r34260 r34844  
    275275      // convert the SkyRectangles to Images for output
    276276      sky_subdivide_image (&image[j*NX_SUB*NY_SUB], &ring[j], NX_SUB, NY_SUB);
    277       printf("%s %8.2f %8.2f\n", ring[j].name, ring[j].coords.crval1, ring[j].coords.crval2);
     277      // printf("%s %8.2f %8.2f\n", ring[j].name, ring[j].coords.crval1, ring[j].coords.crval2);
    278278    }
    279279
    … …  
    691691    // fprintf (stdout, "%f %f  %f  %f  %f %f  %f %f  %f %f %f\n", dec, decUpper, dRA, arg, root1, root2, theta1, theta2, r1, r2, k*r2);
    692692  }
    693   fprintf (stdout, "%f %f  %f x %f (%d)\n", dec, decUpper, dRA, theta, nRA);
     693  // fprintf (stdout, "%f %f  %f x %f (%d)\n", dec, decUpper, dRA, theta, nRA);
    694694
    695695  // I think we need to return the value of dec for the next ring, but I am not sure...
    … …  
    776776      ring[N].coords.crval2 = (j == 0) ? dec / d2r : -dec / d2r;
    777777
    778       printf(" \t %d   %25.20f   %25.20f\n", i, ring[N].coords.crval2, ring[N].coords.crval1);
     778      // printf(" \t %d   %25.20f   %25.20f\n", i, ring[N].coords.crval2, ring[N].coords.crval1);
    779779
    780780      ring[N].coords.pc1_1 = +1.0 * X_PARITY;

  • trunk/Ohana/src/delstar/include/delstar.h

    r34749 r34844  
    2727double ALPHA;
    2828int    VERBOSE;
     29int    VERBOSE2;
    2930int    ORPHAN;
    3031int    MISSED;
    … …  
    3839int          HOST_ID;
    3940char        *HOSTDIR;
     41
     42int SKIP_IMAGES;
    4043
    4144time_t    START;

  • trunk/Ohana/src/delstar/src/args.c

    r34749 r34844  
    4444  }
    4545
     46  VERBOSE2 = FALSE;
     47  if ((N = get_argument (argc, argv, "-vv"))) {
     48    VERBOSE = TRUE;
     49    VERBOSE2 = TRUE;
     50    remove_argument (N, &argc, argv);
     51  }
     52
    4653  // region of interest
    4754  UserPatch.Rmin = 0;
    … …  
    5865    remove_argument (N, &argc, argv);
    5966    UserPatch.Dmax = atof (argv[N]);
     67    remove_argument (N, &argc, argv);
     68  }
     69
     70  SKIP_IMAGES = FALSE;
     71  if ((N = get_argument (argc, argv, "-skip-images"))) {
     72    SKIP_IMAGES = TRUE;
    6073    remove_argument (N, &argc, argv);
    6174  }
    … …  
    205218  }
    206219
     220  VERBOSE2 = FALSE;
     221  if ((N = get_argument (argc, argv, "-v"))) {
     222    VERBOSE = TRUE;
     223    VERBOSE2 = TRUE;
     224    remove_argument (N, &argc, argv);
     225  }
     226
    207227  /* specify portion of the sky */
    208228  UserPatch.Rmin = 0;

  • trunk/Ohana/src/delstar/src/delete_photcodes.c

    r34749 r34844  
    4646
    4747    // an error exit status here is a significant error
    48     if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE, "a")) {
     48    if (!dvo_catalog_open (&catalog, skylist[0].regions[i], VERBOSE2, "a")) {
    4949      fprintf (stderr, "ERROR: failure to open catalog file %s\n", catalog.filename);
    5050      exit (2);
    5151    }
    5252    if (!catalog.Naves_disk) {
    53       if (VERBOSE) fprintf (stderr, "no data in %s, skipping\n", catalog.filename);
     53      if (VERBOSE2) fprintf (stderr, "no data in %s, skipping\n", catalog.filename);
    5454      dvo_catalog_unlock (&catalog);
    5555      dvo_catalog_free (&catalog);
    … …  
    5858
    5959    delete_photcodes_catalog (&catalog, photcodes, Nphotcodes);
    60     dvo_catalog_save_complete (&catalog, VERBOSE);
     60    dvo_catalog_save_complete (&catalog, VERBOSE2);
    6161    dvo_catalog_unlock (&catalog);
    6262    dvo_catalog_free (&catalog);
    … …  
    9797    char tmpline[1024];
    9898    if (VERBOSE)       { snprintf (tmpline, 1024, "%s -v",              command);                    strcpy (command, tmpline); }
     99    if (VERBOSE2)      { snprintf (tmpline, 1024, "%s -vv",             command);                    strcpy (command, tmpline); }
    99100
    100101    fprintf (stderr, "command: %s\n", command);
    … …  
    287288  FREE(measure);
    288289  FREE(average);
     290  FREE(secfilt);
    289291  catalog[0].measure = measureOut;
    290292  catalog[0].average = averageOut;
    … …  
    299301  for (i = 0; i < NaveOut; i++) {
    300302    NmeasureTotal += catalog[0].average[i].Nmeasure;
    301     if (VERBOSE && !(NmeasureTotal <= catalog[0].Nmeasure)) {
     303    if (VERBOSE2 && !(NmeasureTotal <= catalog[0].Nmeasure)) {
    302304      fprintf (stderr, "too few measurements: %d %d %d\n", (int) i, NmeasureTotal, (int) catalog[0].Nmeasure);
    303305    }
    304306    measureOffsetOK &= (catalog[0].average[i].measureOffset < catalog[0].Nmeasure);
    305     if (VERBOSE && !(catalog[0].average[i].measureOffset < catalog[0].Nmeasure)) {
     307    if (VERBOSE2 && !(catalog[0].average[i].measureOffset < catalog[0].Nmeasure)) {
    306308      fprintf (stderr, "offset too large: %d %d %d\n", (int) i, catalog[0].average[i].Nmeasure, (int) catalog[0].Nmeasure);
    307309    }
    308310    measureOffsetOK &= (catalog[0].average[i].measureOffset + catalog[0].average[i].Nmeasure <= catalog[0].Nmeasure);
    309     if (VERBOSE && !(catalog[0].average[i].measureOffset + catalog[0].average[i].Nmeasure <= catalog[0].Nmeasure)) {
     311    if (VERBOSE2 && !(catalog[0].average[i].measureOffset + catalog[0].average[i].Nmeasure <= catalog[0].Nmeasure)) {
    310312      fprintf (stderr, "orrset + Nmeasure too large: %d + %d > %d %d\n", (int) i, catalog[0].average[i].measureOffset, catalog[0].average[i].Nmeasure, (int) catalog[0].Nmeasure);
    311313    }
    … …  
    313315    for (j = 0; j < catalog[0].average[i].Nmeasure; j++) {
    314316      averefOK &= (catalog[0].measure[m+j].averef == i);
    315       if (VERBOSE && !(catalog[0].measure[m+j].averef == i)) {
     317      if (VERBOSE2 && !(catalog[0].measure[m+j].averef == i)) {
    316318        fprintf (stderr, "averef broken: %d vs %d (measure %d)\n", (int) i, catalog[0].measure[m+j].averef, (int) (m+j));
    317319      }
    … …  
    340342  FREE (measureSeqRaw);
    341343  FREE (measureSeqOut);
     344  FREE (measureRefOut);
    342345  FREE (measureAveRaw);
    343346  FREE (averageNmeas);

  • trunk/Ohana/src/delstar/src/delstar.c

    r34749 r34844  
    1010  args (argc, argv);
    1111
    12   set_db (&db);
    13   status = dvo_image_lock (&db, ImageCat, 60.0, LCK_XCLD); /* XCLD */
    14   if (!status) Shutdown ("ERROR: failure to lock image catalog %s", db.filename);
    15   if (db.dbstate == LCK_EMPTY) Shutdown ("ERROR: No images in catalog %s (1)", db.filename);
    16 
    17   status = dvo_image_load (&db, VERBOSE, FALSE);
    18   if (!status) Shutdown ("can't read image catalog %s", db.filename);
     12  if (!SKIP_IMAGES) {
     13    set_db (&db);
     14    status = dvo_image_lock (&db, ImageCat, 60.0, LCK_XCLD); /* XCLD */
     15    if (!status) Shutdown ("ERROR: failure to lock image catalog %s", db.filename);
     16    if (db.dbstate == LCK_EMPTY) Shutdown ("ERROR: No images in catalog %s (1)", db.filename);
     17    status = dvo_image_load (&db, VERBOSE, FALSE);
     18    if (!status) Shutdown ("can't read image catalog %s", db.filename);
     19  }
    1920
    2021  switch (MODE) {
    … …  
    2930      break;
    3031    case MODE_PHOTCODES: {
    31       if (!delete_image_photcodes (&db)) exit (1);
     32      if (!SKIP_IMAGES) {
     33        if (!delete_image_photcodes (&db)) exit (1);
     34      }
    3235
    3336      if (!delete_photcodes (0, NULL)) {

  • trunk/Ohana/src/delstar/src/delstar_client.c

    r34749 r34844  
    3232      break;
    3333    case MODE_PHOTCODES:
    34       delete_photcodes (HOST_ID, HOSTDIR);
    35       break;
     34      if (!delete_photcodes (HOST_ID, HOSTDIR)) {
     35        fprintf (stderr, "failure deleting measurements from %s\n", HOSTDIR);
     36        exit (1);
     37      }
     38      exit (0);
    3639    case MODE_ORPHAN:
    3740      break;

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

    r34749 r34844  
    322322  int *NBAND;
    323323
    324   double   **ra;
    325   double  **dec;
    326   int    **cell;
    327   char  ***name;
     324  double   **ra; // RA of projection cell center
     325  double  **dec; // DEC of projection cell center
     326  int    **cell; // zone,band -> proj cell sequence
     327  char  ***name; // projection cell name
     328 
     329  int NX_SUB;
     330  int NY_SUB;
     331  double dPix;
     332
     333  double **Xo;
     334  double **Yo;
     335  int **dX;
     336  int **dY;
    328337} BoundaryTree;
    329338
    … …  
    496505float PhotKronInst (Measure *measure);
    497506float PhotKronAve (PhotCode *code, Average *average, SecFilt *secfilt);
     507float PhotKronAveErr (PhotCode *code, Average *average, SecFilt *secfilt);
    498508
    499509float PhotZeroPoint (Measure *measure, Average *average, SecFilt *secfilt);
    … …  
    701711int BoundaryTreeSave(char *filename, BoundaryTree *tree);
    702712BoundaryTree *BoundaryTreeLoad(char *filename);
     713int BoundaryTreeProjection (double *x, double *y, double r, double d, BoundaryTree *tree, int zone, int band);
    703714
    704715void dvo_average_init (Average *average);

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

    r34291 r34844  
    4747  gfits_scan (&header, "DEC_ORI", "%lf", 1, &tree->DEC_origin);
    4848  gfits_scan (&header, "DEC_OFF", "%lf", 1, &tree->DEC_offset);
     49
     50  gfits_scan (&header, "NX_SUB", "%d", 1, &tree->NX_SUB);
     51  gfits_scan (&header, "NY_SUB", "%d", 1, &tree->NY_SUB);
     52  gfits_scan (&header, "PIXSCALE", "%lf", 1, &tree->dPix);
    4953
    5054  ftable.header = &theader;
    … …  
    7579    ALLOCATE (tree->ra,   double *, tree->Nzone);
    7680    ALLOCATE (tree->dec,  double *, tree->Nzone);
    77     ALLOCATE (tree->cell, int *, tree->Nzone);
    78     ALLOCATE (tree->name, char **, tree->Nzone);
     81    ALLOCATE (tree->Xo,   double *, tree->Nzone);
     82    ALLOCATE (tree->Yo,   double *, tree->Nzone);
     83    ALLOCATE (tree->dX,      int *, tree->Nzone);
     84    ALLOCATE (tree->dY,      int *, tree->Nzone);
     85    ALLOCATE (tree->cell,    int *, tree->Nzone);
     86    ALLOCATE (tree->name,  char **, tree->Nzone);
    7987    for (i = 0; i < tree->Nzone; i++) {
    8088      ALLOCATE (tree->ra[i],   double, tree->Nband[i]);
    8189      ALLOCATE (tree->dec[i],  double, tree->Nband[i]);
    82       ALLOCATE (tree->cell[i], int,    tree->Nband[i]);
     90      ALLOCATE (tree->Xo[i],   double, tree->Nband[i]);
     91      ALLOCATE (tree->Yo[i],   double, tree->Nband[i]);
     92      ALLOCATE (tree->dX[i],      int, tree->Nband[i]);
     93      ALLOCATE (tree->dY[i],      int, tree->Nband[i]);
     94      ALLOCATE (tree->cell[i],    int, tree->Nband[i]);
    8395      ALLOCATE (tree->name[i], char *, tree->Nband[i]);
    8496      for (j = 0; j < tree->Nband[i]; j++) {
    … …  
    102114    GET_COLUMN_NEW(band,  "BAND",        int);
    103115    GET_COLUMN_NEW(index, "INDEX",       int);
     116    GET_COLUMN_NEW(Xo,    "X_CENT",      double);
     117    GET_COLUMN_NEW(Yo,    "Y_CENT",      double);
     118    GET_COLUMN_NEW(dX,    "X_GRID",      int);
     119    GET_COLUMN_NEW(dY,    "Y_GRID",      int);
    104120    GET_COLUMN_NEW(name,  "NAME",        char); // XXX how is this done?
    105121    gfits_free_header (&theader);
    … …  
    113129      tree->ra[nz][nb] = R[i];
    114130      tree->dec[nz][nb] = D[i];
     131      tree->Xo[nz][nb] = Xo[i];
     132      tree->Yo[nz][nb] = Yo[i];
     133      tree->dX[nz][nb] = dX[i];
     134      tree->dY[nz][nb] = dY[i];
    115135      tree->cell[nz][nb] = i; // XXX ?
    116136      memcpy(tree->name[nz][nb], &name[i*BOUNDARY_TREE_NAME_LENGTH], BOUNDARY_TREE_NAME_LENGTH);
    … …  
    121141    free (zone  );
    122142    free (band  );
     143    free (Xo    );
     144    free (Yo    );
     145    free (dX    );
     146    free (dY    );
    123147    free (index );
    124148    free (name  );
    … …  
    165189  gfits_modify (&header, "DEC_ORI", "%lf", 1, tree->DEC_origin);
    166190  gfits_modify (&header, "DEC_OFF", "%lf", 1, tree->DEC_offset);
     191
     192  gfits_modify (&header, "NX_SUB", "%d", 1, tree->NX_SUB);
     193  gfits_modify (&header, "NY_SUB", "%d", 1, tree->NY_SUB);
     194  gfits_modify (&header, "PIXSCALE", "%lf", 1, tree->dPix);
    167195
    168196  gfits_fwrite_header  (f, &header);
    … …  
    223251    gfits_define_bintable_column (&theader, "J", "BAND", "band sequence number", "none", 1.0, 0.0);
    224252    gfits_define_bintable_column (&theader, "J", "INDEX","cell index", "none", 1.0, 0.0);
     253    gfits_define_bintable_column (&theader, "D", "X_CENT", "projection cell center pixel", "none", 1.0, 0.0);
     254    gfits_define_bintable_column (&theader, "D", "Y_CENT", "projection cell center pixel", "none", 1.0, 0.0);
     255    gfits_define_bintable_column (&theader, "J", "X_GRID", "skycell grid spacing", "none", 1.0, 0.0);
     256    gfits_define_bintable_column (&theader, "J", "Y_GRID", "skycell grid spacing", "none", 1.0, 0.0);
    225257    gfits_define_bintable_column (&theader, fmt, "NAME", "cell name", "none", 1.0, 0.0);
    226258
    … …  
    240272    int    *band          ; ALLOCATE (band,  int,    Ncell);
    241273    int    *index         ; ALLOCATE (index, int,    Ncell);
     274    double *Xo            ; ALLOCATE (Xo,    double, Ncell);
     275    double *Yo            ; ALLOCATE (Yo,    double, Ncell);
     276    int    *dX            ; ALLOCATE (dX,    int,    Ncell);
     277    int    *dY            ; ALLOCATE (dY,    int,    Ncell);
    242278    char   *name          ; ALLOCATE (name,  char,   Ncell*BOUNDARY_TREE_NAME_LENGTH);
    243279
    … …  
    251287        R[i]     = tree->ra[nz][nb];
    252288        D[i]     = tree->dec[nz][nb];
     289        Xo[i]    = tree->Xo[nz][nb];
     290        Yo[i]    = tree->Yo[nz][nb];
     291        dX[i]    = tree->dX[nz][nb];
     292        dY[i]    = tree->dY[nz][nb];
    253293        zone[i]  = nz;
    254294        band[i]  = nb;
    … …  
    265305    gfits_set_bintable_column (&theader, &ftable, "BAND",  band,  Ncell);
    266306    gfits_set_bintable_column (&theader, &ftable, "INDEX", index, Ncell);
     307    gfits_set_bintable_column (&theader, &ftable, "X_CENT", Xo,   Ncell);
     308    gfits_set_bintable_column (&theader, &ftable, "Y_CENT", Yo,   Ncell);
     309    gfits_set_bintable_column (&theader, &ftable, "X_GRID", dX,   Ncell);
     310    gfits_set_bintable_column (&theader, &ftable, "Y_GRID", dY,   Ncell);
    267311    gfits_set_bintable_column (&theader, &ftable, "NAME",  name,  Ncell);
    268312
    … …  
    272316    free (band  );
    273317    free (index );
     318    free (Xo    );
     319    free (Yo    );
     320    free (dX    );
     321    free (dY    );
    274322    free (name  );
    275323
    … …  
    318366}
    319367
     368
     369// projection = TAN
     370// need Ro, Do, Xo, Yo, dPix
     371
     372int BoundaryTreeProjection (double *x, double *y, double r, double d, BoundaryTree *tree, int zone, int band) {
     373
     374    double Xo = tree->Xo[zone][band];
     375    double Yo = tree->Yo[zone][band];
     376    double Ro = tree->ra[zone][band];
     377    double Do = tree->dec[zone][band];
     378    double dPix = tree->dPix;
     379
     380    // this block only depends on Ro, Do
     381
     382    double sdp  = sin(RAD_DEG*Do);
     383    double cdp  = cos(RAD_DEG*Do);
     384    double salp = sin(RAD_DEG*(r - Ro));
     385    double calp = cos(RAD_DEG*(r - Ro));
     386    double sdel = sin(RAD_DEG*d);
     387    double cdel = cos(RAD_DEG*d);
     388   
     389    double stht = sdel*sdp + cdel*cdp*calp;    /* sin(theta) */
     390    double sphi = cdel*salp;                   /* = cos(theta)*sin(phi) */
     391    double cphi = cdel*sdp*calp - sdel*cdp;    /* = cos(theta)*cos(phi) */
     392
     393    // defines the TAN projection (one of zenithal projections available, libdvo/src/coordops.c
     394    // R = cot (theta) = cos(theta) / sin(theta)
     395    double L, M;
     396    if (stht == 0) {
     397        double Rc = hypot(sphi, cphi);
     398        L = 180.0 * sphi / Rc;
     399        M = 180.0 * cphi / Rc;
     400    } else {
     401        L = +DEG_RAD * sphi / stht;
     402        M = -DEG_RAD * cphi / stht;
     403    }
     404
     405    // scale, rotation, parity:
     406    // rotation == 0.0 (pc1_1 == pc2_2 == 1.0, pc1_2 = pc2_1 = 0.0)
     407
     408    // if there were rotation or parity:
     409    // Ro = (coords[0].pc1_1*coords[0].pc2_2 - coords[0].pc1_2*coords[0].pc2_1);
     410    // Xo = (coords[0].pc2_2*L - coords[0].pc1_2*M) / Ro;
     411    // Yo = (coords[0].pc1_1*M - coords[0].pc2_1*L) / Ro;
     412
     413    double X = L;
     414    double Y = M;
     415
     416    // scale is dPix
     417
     418    *x = Xo - X / dPix;
     419    *y = Yo + Y / dPix;
     420   
     421    return TRUE;
     422}
     423

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

    r34620 r34844  
    233233          break;
    234234
     235        case MAG_KRON_ERR:
     236          value.Flt = PhotKronAveErr (field->photcode, average, secfilt);
     237          break;
     238
    235239        case MAG_20:
    236240          value.Flt = PhotM20 (field->photcode, average, secfilt);

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

    r34749 r34844  
    863863    first  = 0;                    // first row in memory to write
    864864    start  = catalog[0].Nsecf_off; // first disk row to write
    865     Nitems = catalog[0].Nsecf_mem;
     865    Nitems = catalog[0].Naverage*Nsecfilt;
    866866    Nrows  = Nitems - first;
    867867

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

    r34620 r34844  
    572572}
    573573
     574float PhotKronAveErr (PhotCode *code, Average *average, SecFilt *secfilt) {
     575
     576  int Ns;
     577  float dMkron;
     578
     579  if (code == NULL) return NAN;
     580
     581  Ns = photcodes[0].hashNsec[code[0].code];
     582  dMkron = (Ns == -1) ? NAN : secfilt[Ns].dMkron;
     583  return (dMkron);
     584}
     585
    574586float PhotAveFluxPSF (PhotCode *code, Average *average, SecFilt *secfilt) {
    575587

  • trunk/Ohana/src/libohana/src

  • trunk/Ohana/src/opihi

  • trunk/Ohana/src/opihi/cmd.astro

  • trunk/Ohana/src/opihi/cmd.data

  • trunk/Ohana/src/opihi/dvo/gimages.c

    r34260 r34844  
    212212        status = XY_to_RD (&Rout, &Dout, X, Y, &image[i].coords);
    213213        fprintf (stderr, "r,d = %f,%f\n", Rout, Dout);
     214        status = RD_to_XY (&X, &Y, Ra, Dec, &image[i].coords);
     215        fprintf (stderr, "x,y = %f,%f\n", X, Y);
    214216    }
    215217# endif

  • trunk/Ohana/src/opihi/dvo/gstar.c

    r34405 r34844  
    11# include "dvoshell.h"
     2
     3typedef enum {
     4  GSTAR_AVE_MAG,
     5  GSTAR_AVE_MAG_ERR,
     6  GSTAR_AVE_MAG_CHISQ,
     7  GSTAR_FILTER_NAMES,
     8  GSTAR_SECF_FLAGS,
     9  GSTAR_AVE_AP_MAG,
     10  GSTAR_AVE_MAG_20,
     11  GSTAR_AVE_MAG_80,
     12  GSTAR_AVE_KRON_MAG,
     13  GSTAR_AVE_KRON_MAG_ERR,
     14  GSTAR_STACK_FLUX_PSF,
     15  GSTAR_STACK_FLUX_PSF_ERR,
     16  GSTAR_STACK_FLUX_KRON,
     17  GSTAR_STACK_FLUX_KRON_ERR,
     18} GSTAR_SECF_CODES;
    219
    320void initPhotcodeSequence (int Nsecfilt);
    … …  
    1330  off_t i, Nstars, *N1;
    1431  off_t j, k, m, N, Nlo, Nhi;
    15   int Nsecfilt, NPTS, QUIET, FULL_OUTPUT, INST;
     32  int Nsecfilt, NPTS, QUIET, FULL_OUTPUT, STACK_OUTPUT, INST;
    1633  int found, GetMeasures;
    1734  int SaveVectors;
    … …  
    5572  if ((N = get_argument (argc, argv, "-full"))) {
    5673    FULL_OUTPUT = TRUE;
     74    remove_argument (N, &argc, argv);
     75  }
     76
     77  STACK_OUTPUT = FALSE;
     78  if ((N = get_argument (argc, argv, "-stack"))) {
     79    STACK_OUTPUT = TRUE;
    5780    remove_argument (N, &argc, argv);
    5881  }
    … …  
    216239     
    217240        /* filter names */
    218         for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, 3);
     241        for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_FILTER_NAMES);
    219242        gprint (GP_LOG, "\n");
    220243
    221244        /* average mags */
    222         for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, 0);
     245        for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_AVE_MAG);
    223246        gprint (GP_LOG, "\n");
    224247
    225248        /* average mag errors */
    226         for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, 1);
     249        for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_AVE_MAG_ERR);
    227250        gprint (GP_LOG, "\n");
    228251
    229252        /* average mag chisq */
    230         for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, 2);
     253        for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_AVE_MAG_CHISQ);
    231254        gprint (GP_LOG, "\n");
    232255
    233256        if (FULL_OUTPUT) {
    234257          /* Map */
    235           for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, 5);
     258          for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_AVE_AP_MAG);
     259          gprint (GP_LOG, "\n");
     260
     261          /* Mkron */
     262          for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_AVE_KRON_MAG);
     263          gprint (GP_LOG, "\n");
     264
     265          /* dMkron */
     266          for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_AVE_KRON_MAG_ERR);
    236267          gprint (GP_LOG, "\n");
    237268
    238269          /* M_20 */
    239           for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, 6);
     270          for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_AVE_MAG_20);
    240271          gprint (GP_LOG, "\n");
    241272
    242273          /* M_80 */
    243           for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, 7);
     274          for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_AVE_MAG_80);;
    244275          gprint (GP_LOG, "\n");
    245276
    246277          /* secfilt flags */
    247           for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, 4);
     278          for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_SECF_FLAGS);
     279          gprint (GP_LOG, "\n");
     280
     281        }
     282        if (STACK_OUTPUT) {
     283          /* FluxPSF */
     284          for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_STACK_FLUX_PSF);
     285          gprint (GP_LOG, "\n");
     286
     287          /* dFluxPSF */
     288          for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_STACK_FLUX_PSF_ERR);
     289          gprint (GP_LOG, "\n");
     290
     291          /* FluxKron */
     292          for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_STACK_FLUX_KRON);
     293          gprint (GP_LOG, "\n");
     294
     295          /* dFluxKron */
     296          for (j = 0; j < Nsecfilt; j++) printPhotcodeSequence (&catalog.average[k], &catalog.secfilt[Nsecfilt*k], j, GSTAR_STACK_FLUX_KRON_ERR);
    248297          gprint (GP_LOG, "\n");
    249298
    … …  
    279328                gprint (GP_LOG, "%6.3f ", catalog.measure[m].Mcal);
    280329                gprint (GP_LOG, "%6.3f ", catalog.measure[m].Map);
     330                gprint (GP_LOG, "%6.3f ", catalog.measure[m].Mkron);
     331                gprint (GP_LOG, "%6.3f ", catalog.measure[m].dMkron);
    281332                gprint (GP_LOG, "%5.1f ", pow(10.0, 0.4*catalog.measure[m].dt));
    282333                gprint (GP_LOG, "%5.3f ", catalog.measure[m].airmass);
    … …  
    296347                gprint (GP_LOG, "%.1f ", catalog.measure[m].extNsigma);
    297348                gprint (GP_LOG, "%.1f ", FromShortDegrees(catalog.measure[m].theta));
     349            }
     350            if (STACK_OUTPUT) {
     351                gprint (GP_LOG, "%10.3e ", catalog.measure[m].FluxPSF);
     352                gprint (GP_LOG, "%10.3e ", catalog.measure[m].dFluxPSF);
     353                gprint (GP_LOG, "%10.3e ", catalog.measure[m].FluxKron);
     354                gprint (GP_LOG, "%10.3e ", catalog.measure[m].dFluxKron);
    298355            }
    299356            gprint (GP_LOG, "\n");
    … …  
    351408}
    352409
     410void print_double_exp (double value) {
     411  if (isnan(value))
     412    gprint (GP_LOG, "  NaN    ");
     413  else
     414    gprint (GP_LOG, "%10.3e ", value);
     415}
     416
    353417void print_short (double value, short int ival) {
    354418  if (ival == NAN_S_SHORT)
    … …  
    386450}
    387451
    388 
    389452void printPhotcodeSequence (Average *average, SecFilt *secfilt, int entry, int type) {
    390453
    … …  
    395458
    396459  switch (type) {
    397     case 0: /* average mags */
     460    case GSTAR_AVE_MAG: /* average mags */
    398461      if (seq == -1) {
    399462        print_double (NAN);
    … …  
    403466      break;
    404467
    405     case 1: /* average mags errors */
     468    case GSTAR_AVE_MAG_ERR: /* average mags errors */
    406469      if (seq == -1) {
    407470        print_double (NAN);
    … …  
    411474      break;
    412475
    413     case 2: /* average mag chisq */
     476    case GSTAR_AVE_MAG_CHISQ: /* average mag chisq */
    414477      if (seq == -1) {
    415478        print_short (NAN_S_SHORT, NAN_S_SHORT);
    … …  
    419482      break;
    420483
    421     case 3: /* filter names */
     484    case GSTAR_FILTER_NAMES: /* filter names */
    422485      code = GetPhotcodebyNsec (seq);
    423486      gprint (GP_LOG, "  %-6s ", code[0].name);
    424487      break;
    425488
    426     case 4: /* secfilt flags */
     489    case GSTAR_SECF_FLAGS: /* secfilt flags */
    427490      if (seq == -1) {
    428491        gprint (GP_LOG, "%08x ", 0);
    … …  
    432495      break;
    433496
    434     case 5: /* average ap mags */
     497    case GSTAR_AVE_AP_MAG: /* average ap mags */
    435498      if (seq == -1) {
    436499        print_double (NAN);
    … …  
    440503      break;
    441504
    442     case 6: /* average ap mags */
     505    case GSTAR_AVE_MAG_20: /* average ap mags */
    443506      if (seq == -1) {
    444507        print_double (NAN);
    … …  
    448511      break;
    449512
    450     case 7: /* average ap mags */
    451       if (seq == -1) {
    452         print_double (NAN);
    453       } else {
    454         print_short (0.001*secfilt[seq].M_20, secfilt[seq].M_80);
    455       }
    456       break;
    457 
    458   }
    459 }
     513    case GSTAR_AVE_MAG_80: /* average ap mags */
     514      if (seq == -1) {
     515        print_double (NAN);
     516      } else {
     517        print_short (0.001*secfilt[seq].M_80, secfilt[seq].M_80);
     518      }
     519      break;
     520
     521    case GSTAR_AVE_KRON_MAG: /* average ap mags */
     522      if (seq == -1) {
     523        print_double (NAN);
     524      } else {
     525        print_double (secfilt[seq].Mkron);
     526      }
     527      break;
     528
     529    case GSTAR_AVE_KRON_MAG_ERR: /* average ap mags */
     530      if (seq == -1) {
     531        print_double (NAN);
     532      } else {
     533        print_double (secfilt[seq].dMkron);
     534      }
     535      break;
     536
     537    case GSTAR_STACK_FLUX_PSF: /* average ap mags */
     538      if (seq == -1) {
     539        print_double (NAN);
     540      } else {
     541        print_double_exp (secfilt[seq].FluxPSF);
     542      }
     543      break;
     544
     545    case GSTAR_STACK_FLUX_PSF_ERR: /* average ap mags */
     546      if (seq == -1) {
     547        print_double (NAN);
     548      } else {
     549        print_double_exp (secfilt[seq].dFluxPSF);
     550      }
     551      break;
     552
     553    case GSTAR_STACK_FLUX_KRON: /* average ap mags */
     554      if (seq == -1) {
     555        print_double (NAN);
     556      } else {
     557        print_double_exp (secfilt[seq].FluxKron);
     558      }
     559      break;
     560
     561    case GSTAR_STACK_FLUX_KRON_ERR: /* average ap mags */
     562      if (seq == -1) {
     563        print_double (NAN);
     564      } else {
     565        print_double_exp (secfilt[seq].dFluxKron);
     566      }
     567      break;
     568
     569  }
     570}

  • trunk/Ohana/src/opihi/pcontrol/HostOps.c

    r26411 r34844  
    232232  FreeIOBuffer (&host[0].comms_buffer);
    233233  FREE (host[0].hostname);
    234   FREE (host[0].job);
    235234  FREE (host);
    236235}

  • trunk/Ohana/src/relphot/Makefile

    r34260 r34844  
    4949$(SRC)/setExclusions.$(ARCH).o   \
    5050$(SRC)/setMrelFinal.$(ARCH).o    \
     51$(SRC)/setMrelCatalog.$(ARCH).o          \
    5152$(SRC)/BoundaryTreeOps.$(ARCH).o         \
    5253$(SRC)/write_coords.$(ARCH).o
    … …  
    7879$(SRC)/setExclusions.$(ARCH).o   \
    7980$(SRC)/setMrelFinal.$(ARCH).o    \
     81$(SRC)/setMrelCatalog.$(ARCH).o    \
    8082$(SRC)/BoundaryTreeOps.$(ARCH).o         \
    8183$(SRC)/write_coords.$(ARCH).o

  • trunk/Ohana/src/relphot/include/relphot.h

    r34642 r34844  
    6363
    6464typedef struct {
     65  int Nfew;
     66  int Ncode;
     67  int Nsys;
     68  int Nbad;
     69  int Ncal;
     70  int Nmos;
     71  int Ngrid;
     72
     73  // NOTE: the following arrays are (possibly) pre-allocated and carried down to each
     74  // thread.  The first 3 (list, dlist, wlist) are used in all relphot analyses; the
     75  // others are only used on the final output steps.
     76
     77  double *list;               // list of measure.mag values for a given star
     78  double *dlist;              // mag errors for a star
     79  double *wlist;              // weights to use for mean mags
     80  double *aplist;             // ap mags for a star
     81  double *kronlist;           // kron mags for a star
     82  double *dkronlist;          // kron mag errors
     83
     84  double *psfqf_list;         // psfqf for all filters
     85  double *psfqfperf_list;     // psfqfperf for all filters
     86  double *stargal_list;       // stargal for all filters
     87
     88  StatType stats;
     89  StatType apstats;
     90  StatType kronstats;
     91} SetMrelInfo;
     92
     93typedef struct {
    6594  AverageTiny *average;       // array of (minimal) average data
    6695  MeasureTiny *measure;       // array of (minimal) measure data
    … …  
    113142
    114143char        *BOUNDARY_TREE;
     144
     145int SET_MREL_VERSION;
    115146
    116147double MAG_LIM;
    … …  
    239270Coords       *getCoords           PROTO((off_t meas, int cat));
    240271off_t         getImageEntry       PROTO((off_t meas, int cat));
     272
     273float         getMcal_alt         PROTO((off_t meas, int cat, FlatCorrectionTable *flatcorr, float Xccd, float Yccd));
    241274float         getMcal             PROTO((off_t meas, int cat, FlatCorrectionTable *flatcorr, Catalog *catalog));
    242275float         getMflat            PROTO((off_t meas, int cat, FlatCorrectionTable *flatcorr, Catalog *catalog));
    … …  
    348381int load_tree (char *treefile);
    349382int BoundaryTreePrimaryCell (char *primaryCellName, double ra, double dec);
     383
     384int print_measure_set_alt (Average *average, SecFilt *secfilt, Measure *measure);
     385int setMrel_catalog_alt (Catalog *catalog, int Nc, int pass, FlatCorrectionTable *flatcorr, SetMrelInfo *results, int Nsecfilt);
     386int setMrelAverageExposure (off_t meas, int cat, int pass, FlatCorrectionTable *flatcorr, SetMrelInfo *results, Average *average, AverageTiny *averageT, SecFilt *secfilt, Measure *measure, MeasureTiny *measureT, off_t *found);
     387int setMrelAverageStack (off_t meas, int cat, FlatCorrectionTable *flatcorr, SetMrelInfo *results, Average *average, SecFilt *secfilt, Measure *measure, char *primaryCell);
     388int setGlobalObjStats (Average *average, Measure *measure);
     389
     390void SetMrelInfoInit (SetMrelInfo *results, int allocLists);

  • trunk/Ohana/src/relphot/src/BoundaryTreeOps.c

    r34260 r34844  
    2121  }
    2222 
     23# define USE_PROJECTION_CELL 0
     24# if (USE_PROJECTION_CELL)
    2325  snprintf (primaryCellName, DVO_MAX_PATH, "RINGS.V3.%s", tree->name[zone][band]);
     26# else
     27
     28  // I have ra, dec, and the primary projection cell.  In order to choose the primary skycell,
     29  // I just need to project to ra,dec to X,Y based on the center of the cell and then get the subdivision right.
     30   
     31  double x = 0.0;
     32  double y = 0.0;
     33  BoundaryTreeProjection (&x, &y, ra, dec, tree, zone, band);
     34 
     35  int xi = x / tree->dX[zone][band];
     36  int yi = y / tree->dY[zone][band];
     37  int N = xi + tree->NX_SUB * yi;
     38 
     39  // XXX short-circuit this for now (we should use this code if we make tree have more variable NX,NY values
     40  // char format[24], skycellname[128];
     41  // int Ndigit = (int)(log10(tree->NX_SUB*tree->NY_SUB)) + 1 ;
     42  // snprintf (format, 24, "%s.%%0%dd", tree->name[zone][band], Ndigit);
     43  // snprintf (skycellname, 128, format, N);
     44
     45  snprintf (primaryCellName, DVO_MAX_PATH, "RINGS.V3.%s.%03d", tree->name[zone][band], N);
     46# endif
     47
    2448  return TRUE;
    2549}

  • trunk/Ohana/src/relphot/src/ImageOps.c

    r34260 r34844  
    352352
    353353// returns image.Mcal - ff(x,y)
     354float getMcal_alt (off_t meas, int cat, FlatCorrectionTable *flatcorr, float Xccd, float Yccd) {
     355
     356  off_t i;
     357  float value, offset;
     358
     359  i = MeasureToImage[cat][meas];
     360  if (i == -1) return (NAN);
     361
     362  if (image[i].flags & IMAGE_BAD) return (NAN); 
     363  value = image[i].Mcal;
     364  offset = 0.0;
     365
     366  // to do this, I need to pass in the catalog and flatcorr pointers
     367  int flat_id = image[i].photom_map_id;
     368  if (flat_id) {
     369    offset = FlatCorrectionOffset (flatcorr, flat_id, Xccd, Yccd);
     370  }
     371  value -= offset;
     372
     373  return (value);
     374}
     375
     376// returns image.Mcal - ff(x,y)
    354377short getUbercalDist (off_t meas, int cat) {
    355378
    … …  
    400423  if (i == -1) return FALSE;
    401424
    402   // this is a bit crude: stack image names are of the form
    403   // RINGS.V3.skycell.1495.027.sky.191211.stk.988232.cmf
    404   // the primaryCell has a name of the form RINGS.V3.skycell.1495
     425  // this is a bit crude: stack image names are of the form:
     426  // RINGS.V3.skycell.1495.027.sky.191211.stk.988232.cmf.
     427
     428  // the primaryCell has a name of the form RINGS.V3.skycell.1495 or RINGS.V3.skycell.1495.027
     429  // (if we use projection or skycell as the primary)
    405430
    406431  if (!strncmp(image[i].name, name, strlen(name))) return TRUE;

  • trunk/Ohana/src/relphot/src/StarOps.c

    r34749 r34844  
    22
    33static int Nmax;
    4 
    5 typedef struct {
    6   int Nfew;
    7   int Ncode;
    8   int Nsys;
    9   int Nbad;
    10   int Ncal;
    11   int Nmos;
    12   int Ngrid;
    13 
    14   // NOTE: the following arrays are (possibly) pre-allocated and carried down to each
    15   // thread.  The first 3 (list, dlist, wlist) are used in all relphot analyses; the
    16   // others are only used on the final output steps.
    17 
    18   double *list;               // list of measure.mag values for a given star
    19   double *dlist;              // mag errors for a star
    20   double *wlist;              // weights to use for mean mags
    21   double *aplist;             // ap mags for a star
    22   double *kronlist;           // kron mags for a star
    23   double *dkronlist;          // kron mag errors
    24 
    25   double *psfqf_list;         // psfqf for all filters
    26   double *psfqfperf_list;     // psfqfperf for all filters
    27   double *stargal_list;       // stargal for all filters
    28 } SetMrelInfo;
    294
    305enum {THREAD_RUN, THREAD_DONE};
    … …  
    154129  for (i = 0; i < Ncatalog; i++) {
    155130    // pass == -1 for anything other than the final pass
    156     setMrel_catalog (catalog, i, -1, flatcorr, &results, Nsecfilt);
     131    switch (SET_MREL_VERSION) {
     132      case 0:
     133        setMrel_catalog (catalog, i, -1, flatcorr, &results, Nsecfilt);
     134        break;
     135      case 1:
     136        setMrel_catalog_alt (catalog, i, -1, flatcorr, &results, Nsecfilt);
     137        break;
     138      default:
     139        fprintf (stderr, "invalid setMrel version %d (use 0 or 1)\n", SET_MREL_VERSION);
     140        exit (2);
     141    }
    157142    SetMrelInfoAccum (&summary, &results);
    158143  }
    … …  
    182167
    183168  for (i = 0; i < Ncatalog; i++) {
    184     setMrel_catalog  (catalog, i, pass, flatcorr, &results, Nsecfilt); // XXX add arguments as needed for options
     169    switch (SET_MREL_VERSION) {
     170      case 0:
     171        setMrel_catalog (catalog, i, pass, flatcorr, &results, Nsecfilt);
     172        break;
     173      case 1:
     174        setMrel_catalog_alt (catalog, i, pass, flatcorr, &results, Nsecfilt);
     175        break;
     176      default:
     177        fprintf (stderr, "invalid setMrel version %d (use 0 or 1)\n", SET_MREL_VERSION);
     178        exit (2);
     179    }
    185180    SetMrelInfoAccum (&summary, &results);
    186181  }
    … …  
    286281
    287282    // pass == -1 for anything other than the final pass
    288     setMrel_catalog (catalog, i, -1, flatcorr, &results, Nsecfilt);
     283    switch (SET_MREL_VERSION) {
     284      case 0:
     285        setMrel_catalog (catalog, i, -1, flatcorr, &results, Nsecfilt);
     286        break;
     287      case 1:
     288        setMrel_catalog_alt (catalog, i, -1, flatcorr, &results, Nsecfilt);
     289        break;
     290      default:
     291        fprintf (stderr, "invalid setMrel version %d (use 0 or 1)\n", SET_MREL_VERSION);
     292        exit (2);
     293    }
    289294    SetMrelInfoAccum (&threadinfo->summary, &results);
    290295  }
    … …  
    332337  liststats_setmode (&kronstats, STATMODE);
    333338
    334   // StatType psfstats;
    335   // liststats_setmode (&psfstats, "MEAN");
    336 
    337339  double *list      = results->list;
    338340  double *dlist     = results->dlist;
    … …  
    363365
    364366    if (isSetMrelFinal) {
     367      // set the name of the primary skycell (this is used in a strcmp to match the skycells in stack detections)
    365368      BoundaryTreePrimaryCell(primaryCell, catalog[Nc].average[j].R, catalog[Nc].average[j].D);
    366369    }

  • trunk/Ohana/src/relphot/src/args.c

    r34642 r34844  
    6464  }
    6565
     66  SET_MREL_VERSION = 1;
     67  if ((N = get_argument (argc, argv, "-set-mrel-version"))) {
     68    remove_argument (N, &argc, argv);
     69    SET_MREL_VERSION = atof(argv[N]);
     70    remove_argument (N, &argc, argv);
     71  }
     72
    6673  VERBOSE = VERBOSE2 = FALSE;
    6774  if ((N = get_argument (argc, argv, "-v"))) {

  • trunk/Ohana/src/uniphot/src/initialize_setposangle.c

    r34749 r34844  
    66    fprintf (stderr, "    -v : verbose mode\n");
    77    fprintf (stderr, "    -update : actually write results to detections tables\n");
     8    fprintf (stderr, "    -parallel : run in parallel mode\n");
     9    fprintf (stderr, "    -h     : this help list\n");
     10    fprintf (stderr, "    -help  : this help list\n");
     11    fprintf (stderr, "    --h    : this help list\n");
     12    fprintf (stderr, "    --help : this help list\n");
    813    fprintf (stderr, "    Note that the dvo db can be specified by -D CATDIR (directory)\n");
    914    exit (2);

  • trunk/Ohana/src/uniphot/src/update_catalog_setposangle.c

    r34749 r34844  
    2525      Mosaic *mosaic = getMosaicForImage(seq);
    2626      Coords *coords = &image[seq].coords;
     27      Coords *mosaicCoords = mosaic ? &mosaic->coords : NULL;
    2728
    2829      double Xccd = catalog[0].measure[m].Xccd;
    2930      double Yccd = catalog[0].measure[m].Yccd;
    3031     
    31       setposangle_local_astrometry (&posAngle, &pltScale, Xccd, Yccd, &mosaic->coords, coords);
     32      setposangle_local_astrometry (&posAngle, &pltScale, Xccd, Yccd, mosaicCoords, coords);
    3233
    3334      catalog[0].measure[m].posangle = ToShortDegrees(posAngle);
    … …  
    5253  // calculate the astrometry for the coordinate of interest
    5354  XY_to_LM (&Lx, &Mx, x,       y,       coords);
    54   XY_to_LM (&Po, &Qo, Lx,      Mx,      mosaic);
     55  if (mosaic) {
     56    XY_to_LM (&Po, &Qo, Lx,      Mx,      mosaic);
     57  } else {
     58    Po = Lx;
     59    Qo = Mx;
     60  }
    5561
    5662  XY_to_LM (&Lx, &Mx, x + 1.0, y,       coords);
    57   XY_to_LM (&Px, &Qx, Lx,      Mx,      mosaic);
     63  if (mosaic) {
     64    XY_to_LM (&Px, &Qx, Lx,      Mx,      mosaic);
     65  } else {
     66    Px = Lx;
     67    Qx = Mx;
     68  }
    5869
    5970  XY_to_LM (&Lx, &Mx, x,       y + 1.0, coords);
    60   XY_to_LM (&Py, &Qy, Lx,      Mx,      mosaic);
     71  if (mosaic) {
     72    XY_to_LM (&Py, &Qy, Lx,      Mx,      mosaic);
     73  } else {
     74    Py = Lx;
     75    Qy = Mx;
     76  }
    6177
    6278  // XXX units for the resulting Tangent Plane coordinates??
Note: See TracChangeset for help on using the changeset viewer.