IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 12206


Ignore:
Timestamp:
Mar 2, 2007, 4:35:02 PM (19 years ago)
Author:
eugene
Message:

updating coordinate handling (access to intermediate frame L,M)

Location:
branches/dvo-mods-2007-02/Ohana/src
Files:
5 edited

Legend:

Unmodified
Added
Removed

  • branches/dvo-mods-2007-02/Ohana/src/libautocode/def/average.d

    r12012 r12206  
    4646# the DVO object IDs are generated internally and are not equivalent to the PSPS object IDs
    4747# probably need to add position chisq
     48
     49# XXX include the number of measurements used to determine the positional information?

  • branches/dvo-mods-2007-02/Ohana/src/libautocode/def/measure.d

    r12008 r12206  
    4747FIELD imageID_lo,     IMAGE_ID_LO,  unsigned int,   reference to image
    4848
    49 # we may need 2 extra bytes for padding purposes...
    50 FIELD dummy,          dummy,        char[6],        padding
     49# note that with airmass = 1.0 / cos(90 - alt), we have full alt/az representation
     50FIELD az,             AZ,           float,          telescope azimuth
     51
     52# we need extra bytes for padding purposes...
     53FIELD dummy,          dummy,        char[2],        padding

  • branches/dvo-mods-2007-02/Ohana/src/libdvo/include/dvo.h

    r12071 r12206  
    3939    PROJ_GLS, // pseudocyl
    4040    PROJ_PAR, // pseudocyl
    41 } OhanaProjections;
     41} OhanaProjection;
    4242
    4343typedef enum {
    … …  
    4646  PROJ_MODE_ZENITHAL,
    4747  PROJ_MODE_PSEUDOCYL,
    48 } OhanaProjectionModes;
     48} OhanaProjectionMode;
    4949
    5050/* RegImage.flag values */
    … …  
    226226OhanaProjection GetProjection (char *ctype);
    227227int SetProjection (char *ctype, OhanaProjection proj);
    228 OhanaProjectionModes GetProjectionMode (OhanaProjections proj);
     228OhanaProjectionMode GetProjectionMode (OhanaProjection proj);
    229229
    230230char *libdvo_version ();

  • branches/dvo-mods-2007-02/Ohana/src/libdvo/src/coordops.c

    r12069 r12206  
    77}
    88
    9 int XY_to_RD (double *ra, double *dec, double x, double y, Coords *coords) {
    10 
    11   int Zenith1, Zenith2, Zenithal, Polynomial, Cartesian, PseudoCyl;
    12   char *type;
    13   double L, M, X, Y, T, Z, Z2;
    14   double R, sphi, cphi, stht, ctht;
     9int XY_to_LM (double *L, double *M, double x, double y, Coords *coords) {
     10
     11  double X, Y;
     12
     13  /** convert pixel coordinates to cartesian system **/
     14  X = coords[0].cdelt1*(x - coords[0].crpix1);
     15  Y = coords[0].cdelt2*(y - coords[0].crpix2);
     16
     17  *L = (X*coords[0].pc1_1 + Y*coords[0].pc1_2);
     18  *M = (X*coords[0].pc2_1 + Y*coords[0].pc2_2);
     19
     20  /** extra polynomial terms **/
     21  if (coords[0].Npolyterms > 1) {
     22    *L += X*X*coords[0].polyterms[0][0] + X*Y*coords[0].polyterms[1][0] + Y*Y*coords[0].polyterms[2][0];
     23    *M += X*X*coords[0].polyterms[0][1] + X*Y*coords[0].polyterms[1][1] + Y*Y*coords[0].polyterms[2][1];
     24  }
     25  if (coords[0].Npolyterms > 2) {
     26    *L += X*X*X*coords[0].polyterms[3][0] + X*X*Y*coords[0].polyterms[4][0] + X*Y*Y*coords[0].polyterms[5][0] + Y*Y*Y*coords[0].polyterms[6][0];
     27    *M += X*X*X*coords[0].polyterms[3][1] + X*X*Y*coords[0].polyterms[4][1] + X*Y*Y*coords[0].polyterms[5][1] + Y*Y*Y*coords[0].polyterms[6][1];
     28  }
     29
     30  return (TRUE);
     31}
     32
     33int LM_to_RD (double *ra, double *dec, double L, double M, Coords *coords) {
     34
     35  OhanaProjection proj;
     36  OhanaProjectionMode mode;
     37  double R, T, Z, Z2, sphi, cphi, stht, ctht;
    1538  double alpha, delta, salp, calp, sdel, sdp, cdp;
    16  
    17   *ra  = 0;
    18   *dec = 0;
    19   stht = ctht = 1;
    20  
     39
    2140  proj = GetProjection (coords[0].ctype);
    2241  mode = GetProjectionMode (proj);
    … …  
    2443  if (proj == PROJ_MODE_NONE) return (FALSE);
    2544
    26   /** convert pixel coordinates to cartesian system **/
    27   X = coords[0].cdelt1*(x - coords[0].crpix1);
    28   Y = coords[0].cdelt2*(y - coords[0].crpix2);
    29 
    30   L = (X*coords[0].pc1_1 + Y*coords[0].pc1_2);
    31   M = (X*coords[0].pc2_1 + Y*coords[0].pc2_2);
    32 
    33   /** extra polynomial terms **/
    34   if (coords[0].Npolyterms > 1) {
    35     L += X*X*coords[0].polyterms[0][0] + X*Y*coords[0].polyterms[1][0] + Y*Y*coords[0].polyterms[2][0];
    36     M += X*X*coords[0].polyterms[0][1] + X*Y*coords[0].polyterms[1][1] + Y*Y*coords[0].polyterms[2][1];
    37   }
    38   if (coords[0].Npolyterms > 2) {
    39     L += X*X*X*coords[0].polyterms[3][0] + X*X*Y*coords[0].polyterms[4][0] + X*Y*Y*coords[0].polyterms[5][0] + Y*Y*Y*coords[0].polyterms[6][0];
    40     M += X*X*X*coords[0].polyterms[3][1] + X*X*Y*coords[0].polyterms[4][1] + X*Y*Y*coords[0].polyterms[5][1] + Y*Y*Y*coords[0].polyterms[6][1];
    41   }
     45  stht = ctht = 1;
    4246
    4347  /** Locally Cartesian Projections **/
    … …  
    7781        break;
    7882      case PROJ_STG:
    79         stht = (4 - RAD_DEG*R) / (4 + RAD_DEC*R);
     83        stht = (4 - RAD_DEG*R) / (4 + RAD_DEG*R);
    8084        ctht = sqrt (1 - stht*stht);
    8185        break;
    … …  
    9498        ctht = sqrt (1 - stht*stht);
    9599        break;
     100      default:
     101        return (FALSE);
    96102    }
    97103    sdp  = sin(RAD_DEG*coords[0].crval2);
    … …  
    115121 
    116122  /**** Other Conventional Projections ****/
    117   if (mode == PROJ_MODE_PSEUDOCLY) {
     123  if (mode == PROJ_MODE_PSEUDOCYL) {
    118124    switch (proj) {
    119125      case PROJ_AIT:
    … …  
    135141        delta = 3 * DEG_RAD * asin (M/180.0);
    136142        break;
     143      default:
     144        return (FALSE);
    137145    }
    138146    *ra  = alpha + coords[0].crval1;
    … …  
    148156}
    149157
    150 int RD_to_XY (double *x, double *y, double ra, double dec, Coords *coords) {
    151 
    152   char *type;
    153   int i, status, Polynomial, Zenith1, Zenith2, Zenithal, Cartesian, PseudoCyl;
     158int RD_to_LM (double *L, double *M, double ra, double dec, Coords *coords) {
     159
    154160  double phi, theta;
    155   double determ;
    156   double X, Y, L, M, Lo, Mo, dL, dM;
     161  double Lo, Mo;
    157162  double sphi, cphi, stht;
    158163  double salp, calp, sdel, cdel, sdp, cdp;
    159164  double P, A, Rc;
    160 
    161   status = TRUE;
    162   *x = 0;
    163   *y = 0;
    164   type = &coords[0].ctype[4];
    165   L = M = 0;
     165  OhanaProjection proj;
     166  OhanaProjectionMode mode;
     167
     168  *L = *M = 0;
    166169
    167170  proj = GetProjection (coords[0].ctype);
    … …  
    175178      if (mosaic == NULL) return (FALSE);
    176179      RD_to_XY (&Lo, &Mo, ra, dec, mosaic);
    177       L = (Lo - coords[0].crval1);
    178       M = (Mo - coords[0].crval2);
    179     } else {
    180       L = (ra  - coords[0].crval1);
    181       M = (dec - coords[0].crval2);
    182     }
     180      *L = (Lo - coords[0].crval1);
     181      *M = (Mo - coords[0].crval2);
     182      return (TRUE);
     183    }
     184    *L = (ra  - coords[0].crval1);
     185    *M = (dec - coords[0].crval2);
     186    return (TRUE);
    183187  }
    184188 
    … …  
    195199    sphi = cdel*salp;                   /* = cos(theta)*sin(phi) */
    196200    cphi = cdel*sdp*calp - sdel*cdp;    /* = cos(theta)*cos(phi) */
    197     if (stht < 0) status = FALSE;
     201    if (stht < 0) return (FALSE);
    198202
    199203    switch (proj) {
    200204      case PROJ_TAN:
    201205      case PROJ_DIS:
    202         L =  DEG_RAD * sphi / stht;
    203         M = -DEG_RAD * cphi / stht;
    204         break;
     206        *L =  DEG_RAD * sphi / stht;
     207        *M = -DEG_RAD * cphi / stht;
     208        return (TRUE);
    205209      case PROJ_SIN:
    206         L =  DEG_RAD * sphi;
    207         M = -DEG_RAD * cphi;
    208         break;
     210        *L =  DEG_RAD * sphi;
     211        *M = -DEG_RAD * cphi;
     212        return (TRUE);
    209213      case PROJ_ZEA:
    210214      case PROJ_ZPL:
    211215        Rc = DEG_RAD * M_SQRT2 / sqrt (1 + stht);
    212         L =  Rc * sphi;
    213         M = -Rc * cphi;
    214         status = TRUE;
    215         break;
    216     }
     216        *L =  Rc * sphi;
     217        *M = -Rc * cphi;
     218        return (TRUE);
     219      default:
     220        return (FALSE);
     221    }
     222    return (FALSE);
    217223  }
    218224
    … …  
    224230        theta = RAD_DEG*(dec - coords[0].crval2);
    225231        P = 1.0 + cos (theta) * cos (0.5*phi);
    226         if (P != 0.0) {
    227           A =  DEG_RAD * sqrt (2.0 / P);
    228           L =  2.0 * A * cos (theta) * sin (0.5*phi);
    229           M =  A * sin (theta);
    230         } else {
    231           L =  0.0;
    232           M =  0.0;
    233         }       
    234         break;
     232        if (P == 0.0) {
     233          *L =  0.0;
     234          *M =  0.0;
     235          return (TRUE);
     236        }
     237        A =  DEG_RAD * sqrt (2.0 / P);
     238        *L =  2.0 * A * cos (theta) * sin (0.5*phi);
     239        *M =  A * sin (theta);
     240        return (TRUE);
    235241      case PROJ_GLS:
    236242        phi = ra - coords[0].crval1;
    237243        theta = dec - coords[0].crval2;
    238         L = phi * cos(RAD_DEG * theta);
    239         M = theta;
    240         break;
     244        *L = phi * cos(RAD_DEG * theta);
     245        *M = theta;
     246        return (TRUE);
    241247      case PROJ_PAR:
    242248        phi = ra - coords[0].crval1;
    243249        theta = dec - coords[0].crval2;
    244         L = phi * (2.0*cos(2*RAD_DEG*theta/3.0) - 1);
    245         M = 180.0 * sin (RAD_DEG*theta/3.0);
    246         break;
    247     }
    248   }
     250        *L = phi * (2.0*cos(2*RAD_DEG*theta/3.0) - 1);
     251        *M = 180.0 * sin (RAD_DEG*theta/3.0);
     252        return (TRUE);
     253      default:
     254        return (FALSE);
     255    }
     256    return (FALSE);
     257  }
     258  return (FALSE);
     259}
     260
     261int LM_to_XY (double *x, double *y, double L, double M, Coords *coords) {
     262
     263  int i;
     264  double determ;
     265  double X, Y, Lo, Mo, dL, dM;
     266
     267  *x = 0;
     268  *y = 0;
    249269
    250270  /* convert L,M to X,Y */
    … …  
    278298  *y = Y / coords[0].cdelt2 + coords[0].crpix2;
    279299
     300  return (TRUE);
     301}
     302
     303int XY_to_RD (double *ra, double *dec, double x, double y, Coords *coords) {
     304
     305  double L, M;
     306  int status;
     307
     308  status = XY_to_LM (&L, &M, x, y, coords);
     309  if (!status) return FALSE;
     310
     311  status = LM_to_RD (ra, dec, L, M, coords);
     312  return (status);
     313}
     314
     315int RD_to_XY (double *x, double *y, double ra, double dec, Coords *coords) {
     316
     317  double L, M;
     318  int status;
     319
     320  status = RD_to_LM (&L, &M, ra, dec, coords);
     321  if (!status) return FALSE;
     322
     323  status = LM_to_XY (x, y, L, M, coords);
    280324  return (status);
    281325}
    … …  
    601645     WRP is equiv to PLY, with implied mosaic */
    602646
    603 int GetProjection (char *ctype) {
     647OhanaProjection GetProjection (char *ctype) {
    604648  if (!strcmp(&ctype[4], "-ZEA")) return PROJ_ZEA;
    605649  if (!strcmp(&ctype[4], "-ZPL")) return PROJ_ZPL;
    … …  
    620664}
    621665 
    622 int SetProjection (char *ctype, int proj) {
     666int SetProjection (char *ctype, OhanaProjection proj) {
    623667  switch (proj) {
    624     case PROJ_ZEA: strcpy(&ctype[4], "-ZEA") return TRUE;
    625     case PROJ_ZPL: strcpy(&ctype[4], "-ZPL") return TRUE;
    626     case PROJ_ARC: strcpy(&ctype[4], "-ARC") return TRUE;
    627     case PROJ_STG: strcpy(&ctype[4], "-STG") return TRUE;
    628     case PROJ_SIN: strcpy(&ctype[4], "-SIN") return TRUE;
    629     case PROJ_TAN: strcpy(&ctype[4], "-TAN") return TRUE;
    630     case PROJ_DIS: strcpy(&ctype[4], "-DIS") return TRUE;
    631     case PROJ_LIN: strcpy(&ctype[4], "-LIN") return TRUE;
    632     case PROJ_PLY: strcpy(&ctype[4], "-PLY") return TRUE;
    633     case PROJ_WRP: strcpy(&ctype[4], "-WRP") return TRUE;
    634     case PROJ_AIT: strcpy(&ctype[4], "-AIT") return TRUE;
    635     case PROJ_GLS: strcpy(&ctype[4], "-GLS") return TRUE;
    636     case PROJ_PAR: strcpy(&ctype[4], "-PAR") return TRUE;
     668    case PROJ_ZEA: strcpy(&ctype[4], "-ZEA"); return TRUE;
     669    case PROJ_ZPL: strcpy(&ctype[4], "-ZPL"); return TRUE;
     670    case PROJ_ARC: strcpy(&ctype[4], "-ARC"); return TRUE;
     671    case PROJ_STG: strcpy(&ctype[4], "-STG"); return TRUE;
     672    case PROJ_SIN: strcpy(&ctype[4], "-SIN"); return TRUE;
     673    case PROJ_TAN: strcpy(&ctype[4], "-TAN"); return TRUE;
     674    case PROJ_DIS: strcpy(&ctype[4], "-DIS"); return TRUE;
     675    case PROJ_LIN: strcpy(&ctype[4], "-LIN"); return TRUE;
     676    case PROJ_PLY: strcpy(&ctype[4], "-PLY"); return TRUE;
     677    case PROJ_WRP: strcpy(&ctype[4], "-WRP"); return TRUE;
     678    case PROJ_AIT: strcpy(&ctype[4], "-AIT"); return TRUE;
     679    case PROJ_GLS: strcpy(&ctype[4], "-GLS"); return TRUE;
     680    case PROJ_PAR: strcpy(&ctype[4], "-PAR"); return TRUE;
     681    case PROJ_NONE: return FALSE;
    637682  }
    638683  return FALSE;
    639684} 
    640685
    641 int GetProjectionMode (int proj) {
     686OhanaProjectionMode GetProjectionMode (OhanaProjection proj) {
    642687  switch (proj) {
    643688    case PROJ_ZEA:
    … …  
    656701    case PROJ_GLS:
    657702    case PROJ_PAR:
    658       return PROJ_MODE_PSEUDOCLY;
     703      return PROJ_MODE_PSEUDOCYL;
     704    default: PROJ_MODE_NONE;
    659705  }
    660706  return PROJ_MODE_NONE;

  • branches/dvo-mods-2007-02/Ohana/src/relastro/include/relastro.h

    r12050 r12206  
    33# include <kapa.h>
    44# include <signal.h>
     5
     6// structure to hold coordinate fitting terms
     7typedef struct {
     8    int Npts;
     9    int Nterms;
     10    int Npower;
     11    int Norder;
     12    double **xsum;
     13    double **ysum;
     14
     15    double **vector;
     16    double **matrix;
     17} CoordFit;
    518
    619typedef struct {
    … …  
    1225
    1326  double p, dp;
     27
     28  double chisq;
     29  int Nfit;
    1430} PMfit;
    1531
Note: See TracChangeset for help on using the changeset viewer.