IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 3326


Ignore:
Timestamp:
Feb 24, 2005, 4:53:52 PM (21 years ago)
Author:
eugene
Message:

added PutCoords, cleaned up type values

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/Ohana/src/libohana/src/coordops.c

    r3296 r3326  
    1515  type = &coords[0].ctype[4];
    1616 
     17  /* DIS and PLY are identical & apply projection; WRP is cartesian */
    1718  Polynomial = !strcmp(type, "-PLY") || !strcmp(type, "-DIS") || !strcmp(type, "-WRP");
    18   Zenithal   = !strcmp(type, "-PLY") || !strcmp(type, "-TAN") || !strcmp(type, "-SIN") || !strcmp(type, "-ZEA") || !strcmp(&coords[0].ctype[0], "MM");
     19  Zenithal   = !strcmp(type, "-PLY") || !strcmp(type, "-DIS") || !strcmp(type, "-TAN") || !strcmp(type, "-SIN") || !strcmp(type, "-ZEA") || !strcmp(&coords[0].ctype[0], "MM");
    1920  Cartesian  = !strcmp(type, "-LIN") || !strcmp(type, "-WRP") || !strcmp(&coords[0].ctype[0], "GENE");
    2021  PseudoCyl  = !strcmp(type, "-AIT") || !strcmp(type, "-GLS") || !strcmp(type, "-PAR");
     
    5859    }
    5960
    60     if (!strcmp(type, "-PLY") || !strcmp(type, "-TAN") || !strcmp(type, "-DIS")) {
     61    if (!strcmp(type, "-PLY") || !strcmp(type, "-DIS") || !strcmp(type, "-TAN")) {
    6162      if (R == 0) {
    6263        stht = 1.0;
     
    137138
    138139  Polynomial = !strcmp(type, "-PLY") || !strcmp(type, "-DIS") || !strcmp(type, "-WRP");
    139   Zenithal   = !strcmp(type, "-PLY") || !strcmp(type, "-TAN") || !strcmp(type, "-SIN") || !strcmp(type, "-ZEA") || !strcmp(&coords[0].ctype[0], "MM");
     140  Zenithal   = !strcmp(type, "-PLY") || !strcmp(type, "-DIS") || !strcmp(type, "-TAN") || !strcmp(type, "-SIN") || !strcmp(type, "-ZEA") || !strcmp(&coords[0].ctype[0], "MM");
    140141  Cartesian  = !strcmp(type, "-LIN") || !strcmp(type, "-WRP") || !strcmp(&coords[0].ctype[0], "GENE");
    141142  PseudoCyl  = !strcmp(type, "-AIT") || !strcmp(type, "-GLS") || !strcmp(type, "-PAR");
     
    162163    if (stht < 0) status = FALSE;
    163164
    164     if (!strcmp(type, "-PLY") || !strcmp(type, "-TAN") || !strcmp(type, "-DIS")) {
     165    if (!strcmp(type, "-PLY") || !strcmp(type, "-DIS") || !strcmp(type, "-TAN")) {
    165166      L =  DEG_RAD * sphi / stht;
    166167      M = -DEG_RAD * cphi / stht;
     
    269270int GetCoords (Coords *coords, Header *header) {
    270271 
    271   int status, itmp;
    272   double rotate, scale;
     272  int status, itmp, Polynomial;
     273  double Lambda, rotate, scale;
    273274  double equinox;
     275  char *ctype;
    274276 
    275277  rotate = 0.0;
     
    290292      status &= fits_scan (header, "CDELT2", "%f", 1, &coords[0].cdelt2);
    291293      if (fits_scan (header, "CROTA2", "%f", 1, &rotate)) {
     294        Lambda = coords[0].cdelt2 / coords[0].cdelt1;
    292295        coords[0].pc1_1 =  cos(rotate*RAD_DEG);
    293         coords[0].pc1_2 = -sin(rotate*RAD_DEG);
    294         coords[0].pc2_1 =  sin(rotate*RAD_DEG);
     296        coords[0].pc1_2 = -sin(rotate*RAD_DEG) * Lambda;
     297        coords[0].pc2_1 =  sin(rotate*RAD_DEG) / Lambda;
    295298        coords[0].pc2_2 =  cos(rotate*RAD_DEG);
    296299      }
     
    300303        status &= fits_scan (header, "PC002002", "%f", 1, &coords[0].pc2_2);
    301304      }
    302       if (!strcmp (coords[0].ctype, "DEC--PLY")) {
     305
     306      ctype = &coords[0].ctype[4];
     307      Polynomial = !strcmp (ctype, "-PLY") || !strcmp (ctype, "-DIS") || !strcmp (ctype, "-WRP");
     308
     309      if (Polynomial) {
    303310        if (fits_scan (header, "NPLYTERM", "%d", 1, &itmp)) {
    304311          coords[0].Npolyterms = itmp;
     
    345352    }
    346353  } else {
     354    /* some of my thesis data uses this simple linear model - convert on read? */
    347355    if (fits_scan (header, "RA_O", "%lf", 1, &coords[0].crval1)) {
    348356      status  = fits_scan (header, "RA_X", "%f", 1, &coords[0].pc1_1);
     
    375383  }
    376384  return (status);
     385}
     386
     387int PutCoords (Coords *coords, Header *header) {
     388 
     389  int OldAIPS;
     390  char csys[16], ctype[16];
     391  double rotate, Lambda;
     392
     393  /* modifications to the ctype? */
     394  OldAIPS = FALSE;
     395  fits_modify (header, "CTYPE2",   "%s",  1, coords[0].ctype);
     396  if (!strcmp(coords[0].ctype, "MM")) {
     397    fits_modify (header, "CTYPE1",   "%s",  1, "LL");
     398    OldAIPS = TRUE;
     399  } else {
     400    strcpy (csys, "NONE");
     401    if (!strncmp (coords[0].ctype, "DEC-", 4)) strcpy (csys, "RA--");
     402    if (!strncmp (coords[0].ctype, "GLAT", 4)) strcpy (csys, "GLON");
     403    if (!strncmp (coords[0].ctype, "ELAT", 4)) strcpy (csys, "ELON");
     404    if (!strncmp (coords[0].ctype, "HLAT", 4)) strcpy (csys, "HLON");
     405    if (!strncmp (coords[0].ctype, "SLAT", 4)) strcpy (csys, "SLON");
     406    if (!strcmp (csys, "NONE")) return (FALSE);
     407    sprintf (ctype, "%s-%s", csys, &coords[0].ctype[5]);
     408    fits_modify (header, "CTYPE1",   "%s",  1, ctype);
     409  }   
     410
     411  fits_modify (header, "CDELT1",   "%le", 1, coords[0].cdelt1);
     412  fits_modify (header, "CDELT2",   "%le", 1, coords[0].cdelt2);
     413  fits_modify (header, "CRVAL1",   "%lf", 1, coords[0].crval1);
     414  fits_modify (header, "CRVAL2",   "%lf", 1, coords[0].crval2); 
     415  fits_modify (header, "CRPIX1",   "%lf", 1, coords[0].crpix1);
     416  fits_modify (header, "CRPIX2",   "%lf", 1, coords[0].crpix2);
     417
     418  if (OldAIPS) {
     419    Lambda = coords[0].cdelt2 / coords[0].cdelt1;
     420    rotate = DEG_RAD*atan2 (coords[0].pc2_1*Lambda, coords[0].pc1_1);
     421    fits_modify (header, "CROTA1", "%f", 1, rotate);
     422    fits_modify (header, "CROTA2", "%f", 1, rotate);
     423    return (TRUE);
     424  }
     425
     426  fits_modify (header, "PC001001", "%le", 1, coords[0].pc1_1);
     427  fits_modify (header, "PC001002", "%le", 1, coords[0].pc1_2);
     428  fits_modify (header, "PC002001", "%le", 1, coords[0].pc2_1);
     429  fits_modify (header, "PC002002", "%le", 1, coords[0].pc2_2);
     430  fits_modify (header, "NPLYTERM", "%d",  1, coords[0].Npolyterms);
     431
     432  /* RA Terms */
     433  if (coords[0].Npolyterms > 1) {
     434    fits_modify (header, "PCA1X2Y0", "%le", 1, coords[0].polyterms[0][0]);   /* polyterms[0]); */
     435    fits_modify (header, "PCA1X1Y1", "%le", 1, coords[0].polyterms[1][0]);   /* polyterms[1]); */
     436    fits_modify (header, "PCA1X0Y2", "%le", 1, coords[0].polyterms[2][0]);   /* polyterms[2]); */
     437  }
     438  if (coords[0].Npolyterms > 2) {
     439    fits_modify (header, "PCA1X3Y0", "%le", 1, coords[0].polyterms[3][0]);   /* polyterms[3]); */
     440    fits_modify (header, "PCA1X2Y1", "%le", 1, coords[0].polyterms[4][0]);   /* polyterms[4]); */
     441    fits_modify (header, "PCA1X1Y2", "%le", 1, coords[0].polyterms[5][0]);   /* polyterms[5]); */
     442    fits_modify (header, "PCA1X0Y3", "%le", 1, coords[0].polyterms[6][0]);   /* polyterms[6]); */
     443  }
     444
     445  /* Dec Terms */
     446  if (coords[0].Npolyterms > 1) {
     447    fits_modify (header, "PCA2X2Y0", "%le", 1, coords[0].polyterms[0][1]);   /* polyterms[7]); */
     448    fits_modify (header, "PCA2X1Y1", "%le", 1, coords[0].polyterms[1][1]);   /* polyterms[8]); */
     449    fits_modify (header, "PCA2X0Y2", "%le", 1, coords[0].polyterms[2][1]);   /* polyterms[9]); */
     450  }
     451  if (coords[0].Npolyterms > 2) {
     452    fits_modify (header, "PCA2X3Y0", "%le", 1, coords[0].polyterms[3][1]);   /* polyterms[10]); */
     453    fits_modify (header, "PCA2X2Y1", "%le", 1, coords[0].polyterms[4][1]);   /* polyterms[11]); */
     454    fits_modify (header, "PCA2X1Y2", "%le", 1, coords[0].polyterms[5][1]);   /* polyterms[12]); */
     455    fits_modify (header, "PCA2X0Y3", "%le", 1, coords[0].polyterms[6][1]);   /* polyterms[13]); */
     456  }
     457  return (TRUE);
    377458}
    378459
Note: See TracChangeset for help on using the changeset viewer.