IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

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:
6 edited

Legend:

Unmodified
Added
Removed
  • trunk/Ohana

  • 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
Note: See TracChangeset for help on using the changeset viewer.