IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 3401


Ignore:
Timestamp:
Mar 9, 2005, 9:22:34 PM (21 years ago)
Author:
eugene
Message:

added mosaic astrom support (WRP / DIS)

Location:
trunk/Ohana/src/mosastro
Files:
22 edited

Legend:

Unmodified
Added
Removed

  • trunk/Ohana/src/mosastro/Makefile

    r3323 r3401  
    5656$(SRC)/GetScatter.$(ARCH).o \
    5757$(SRC)/clip.$(ARCH).o \
    58 $(SRC)/parse_time.$(ARCH).o
     58$(SRC)/parse_time.$(ARCH).o \
     59$(SRC)/testcoords.$(ARCH).o
    5960
    6061mosastro: $(BIN)/mosastro.$(ARCH)
    … …  
    102103$(STD) : $(INC)/mosastro.h
    103104
    104 INSTALL = mosastro mkstandards mkobs
     105# INSTALL = mosastro mkstandards mkobs
     106INSTALL = mosastro
    105107
    106108# dependancy rules for binary code #########################

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

    r3323 r3401  
    9696char *FOCAL_PLANE;
    9797
     98CatStats *gregions (CatStats *patch, int *nregion);
     99Gradients *GetGradients ();
     100Header *mkheader (int Nx, int Ny, int Nstars, Coords *coords);
     101Header *mkmosaic (int Nx, int Ny, int Nstars, Coords *coords);
    98102StarData *getusno (CatStats *catstats, int *Nstars);
    99103StarData *getgsc (CatStats *catstats, int *NSTARS);
    100104StarData *get2mass (CatStats *catstats, int *NSTARS);
    101105StarData *getptolemy (CatStats *catstats, int *NSTARS);
    102 void GetConfig (char *config, char *field, char *format, int N, void *ptr);
    103106StarData *greference (int *Nrefcat);
    104 
    105 Header *mkheader (int Nx, int Ny, int Nstars, Coords *coords);
    106 
    107107StarData *gcatalog (char *filename, int *Nstars);
    108 
    109 CatStats *gregions (CatStats *patch, int *nregion);
     108StarData *gptolemy (char *filename, int *NSTARS);
     109double GetScatter ();
    110110int find_dec_bands (CatStats *area);
    111111int load_ra_blocks (int Ndec, CatStats *area);
    112112int parse_GSC_line (CatStats *tregion, char *line);
     113int fake_field (double RA, double DEC);
     114int LoadStars (int Nfile, char **file);
     115int project_refcat (StarData *refcat, int Nrefcat);
     116int deproject_stars ();
     117int project_stars ();
     118int deproject_raw ();
     119int project_ref ();
     120int init_field ();
     121int load_field (char *filename);
     122int init_chips ();
     123int load_chips (char *filename);
     124int dump_stars (FILE *f, StarData *stars, int Nstars);
     125int dump_grads (Gradients *grad, char *filename);
     126int dump_match ();
     127int dump_rawstars ();
     128int dump_refcat (StarData *refcat, int Nrefcat);
     129int match (StarData *refcat, int Nrefcat);
     130int mkvector (int n, int m, int norder);
     131int ClipOnFP (double Nsigma);
     132void GetConfig (char *config, char *field, char *format, int N, void *ptr);
    113133void init_regions ();
    114134void add_to_regions (CatStats *area);
    115135void area_of_region (CatStats *region);
    116136void set_catalog (char *catdir);
    117 
    118 StarData *gptolemy (char *filename, int *NSTARS);
    119 int fake_field (double RA, double DEC);
    120 Gradients *GetGradients ();
    121 void FitChips (int Norder);
    122137void FitChipLinear (StarData *raw, StarData *ref, int Nmatch, Coords *coords);
    123138void FitChipResid (StarData *raw, StarData *ref, int Nmatch, Coords *coords);
    124 double GetScatter ();
    125 Header *mkmosaic (int Nx, int Ny, int Nstars, Coords *coords);
     139void ConfigInit (int *argc, char **argv);
     140void print_help ();
     141void args (int *argc, char **argv);
     142void field_stats ();
     143void field_combine ();
     144void fit_apply_grads (Coords *distort, Coords *project, int term);
     145void fit_apply_coords (Coords *coords);
     146void fit_eval ();
     147void fit_add (double x1, double y1, double x2, double y2);
     148void fit_init (int order);
     149void fit_correct_grads (Gradients *in, Gradients *out, int term);
     150void fit_free ();
     151void output (char *ext, char *phu);
     152void FitGradients (Gradients *grad);
     153void uppercase (char *string);
     154e_time parse_time (Header *header);
     155int sortthree (double *X, double *Y, int *Z, int N);
     156StarData *getstone (CatStats *input, int *nstars);
     157void ChipToSky (StarData *stars, int Nstars, Coords *coords);
     158void ChipToFP (StarData *stars, int Nstars, Coords *coords);
     159void FPtoTP (StarData *stars, int Nstars, Coords *coords);
     160void TPtoSky (StarData *stars, int Nstars, Coords *coords);
     161void SkyToTP (StarData *stars, int Nstars, Coords *coords);
     162void TPtoFP (StarData *stars, int Nstars, Coords *coords);
     163void FPtoChip (StarData *stars, int Nstars, Coords *coords);
     164int gaussj (double **a, int n, double **b, int m);
     165void FitChips (int Norder);
     166void FitChip (StarData *raw, StarData *ref, int Nmatch, Coords *coords);
     167int mkpolyterm (int n, int m);
     168void wchip (char *filename, Chip *data);
     169void wstars (char *filename, Stars *stars, int Nstars, Header *header);

  • trunk/Ohana/src/mosastro/src/FitChips.c

    r3323 r3401  
    1919  fit_init (coords[0].Npolyterms);
    2020  for (i = 0; i < Nmatch; i++) {
     21    if (raw[i].mask) continue;
    2122    fit_add (raw[i].X, raw[i].Y, ref[i].L, ref[i].M);
    2223  }

  • trunk/Ohana/src/mosastro/src/FitGradients.c

    r3323 r3401  
    44void FitGradients (Gradients *grad) {
    55
    6   int i, j, Nterm, Norder;
    7   double **Pf, **Qf, chisq;
     6  int i, Norder;
    87  Gradients grfix;
    98

  • trunk/Ohana/src/mosastro/src/GetScatter.c

    r3323 r3401  
    2121    }
    2222  }
    23   dR = sqrt (dR2 / Ntotal);
     23  dR = sqrt (dR2 / Ntotal) * 3600.0 * field.project.cdelt1;
    2424  *Nscatter = Ntotal;
    2525  return (dR);

  • trunk/Ohana/src/mosastro/src/LoadStars.c

    r3323 r3401  
    44int LoadStars (int Nfile, char **file) {
    55
    6   Header *header;
    7   int i, j, NCHIP, itmp, Nbytes, nbytes;
     6  int i, j, itmp, Nbytes, nbytes;
    87  FILE *f;
    98  time_t tsval;

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

    r3323 r3401  
    1212 
    1313  int N;
    14   char line[500];
    1514
    1615  if (get_argument (*argc, argv, "-help") ||

  • trunk/Ohana/src/mosastro/src/chips.c

    r3323 r3401  
    22
    33/* set chip model based on initial header info */
     4
     5/** when the chips and field terms are treated independently, they must be of
     6    ctype PLY.  When we reach a solution, and write out a coupled set of chip
     7    and field terms, then we need to transition to WRP (chip) & DIS (field
     8**/
     9
    410int init_chips () {
    511
    … …  
    1521
    1622    /* bore site center guess */
    17     strcpy (chip[i].map.ctype, "DEC--WRP");
     23    strcpy (chip[i].map.ctype, "DEC--PLY");
    1824    chip[i].map.crpix1 = 0.0;
    1925    chip[i].map.crpix2 = 0.0;
    … …  
    5359  exit (1);
    5460}
     61

  • trunk/Ohana/src/mosastro/src/clip.c

    r3323 r3401  
    33int ClipOnFP (double Nsigma) {
    44
    5   int i, j, Nscatter;
     5  int i, j, Nscatter, Nmask, Nkeep;
    66  double dL, dM, dV;
    77  double sigma;
    … …  
    99
    1010  sigma = GetScatter (&Nscatter);
     11  fprintf (stderr, "scatter (chip) : %f for %d stars\n", sigma, Nscatter);
     12
     13  Nmask = Nkeep = 0;
    1114
    1215  for (i = 0; i < Nchip; i++) {
    1316    raw = chip[i].raw;
    1417    ref = chip[i].ref;
    15     for (j = 0; j < chip[i].Nstars; j++) {
     18    for (j = 0; j < chip[i].Nmatch; j++) {
    1619      dL = raw[j].L - ref[j].L;
    1720      dM = raw[j].M - ref[j].M;
    … …  
    1922      if (dV > Nsigma*sigma) {
    2023        raw[j].mask = TRUE;
     24        Nmask ++;
    2125      } else {
    2226        raw[j].mask = FALSE;   
     27        Nkeep ++;
    2328      }
    2429    }
    2530  }
     31
     32  fprintf (stderr, "Nmask: %d, Nkeep: %d\n", Nmask, Nkeep);
    2633  return (TRUE);
    2734}

  • trunk/Ohana/src/mosastro/src/dump.c

    r3323 r3401  
    6666int dump_refcat (StarData *refcat, int Nrefcat) {
    6767
    68   int i;
    6968  FILE *f;
    7069  f = fopen ("refcat.dat", "w");

  • trunk/Ohana/src/mosastro/src/fakefield.c

    r3323 r3401  
    3838
    3939  /** distort only has power in polyterms **/
    40   strcpy (field.distort.ctype, "DEC--WRP");
     40  strcpy (field.distort.ctype, "DEC--PLY");
    4141  field.distort.crval1 = 0.0;
    4242  field.distort.crval2 = 0.0;

  • trunk/Ohana/src/mosastro/src/field.c

    r3323 r3401  
    22
    33/* determine an initial guess to field parameters from data */
    4 int init_field (int Norder) {
     4int init_field () {
    55
    6   int i, Nterm;
    7   double PS;
     6  int i;
    87
    98  if (FIELD != (char *) NULL) {
    … …  
    3837
    3938  /* bore site center guess */
    40   strcpy (field.distort.ctype, "DEC--WRP");
     39  strcpy (field.distort.ctype, "DEC--PLY");
    4140  field.distort.crval1 = 0.0;
    4241  field.distort.crval2 = 0.0;
    … …  
    8988
    9089  int i;
     90  double cd1, cd2, pc11, pc12, pc21, pc22;
    9191
    92   strcpy (field.project.ctype, "DEC--WRP");
     92  /* combine boresite & distortion parameters: ctype DIS */
     93  strcpy (field.project.ctype, "DEC--DIS");
     94
     95  cd1  = field.project.cdelt1;
     96  cd2  = field.project.cdelt2;
     97  pc11 = field.project.pc1_1;
     98  pc12 = field.project.pc1_2;
     99  pc21 = field.project.pc2_1;
     100  pc22 = field.project.pc2_2;
    93101
    94102  field.project.Npolyterms = field.distort.Npolyterms;
     103
    95104  for (i = 0; i < 7; i++) {
    96     field.project.polyterms[i][0] = field.distort.polyterms[i][0];
    97     field.project.polyterms[i][1] = field.distort.polyterms[i][1];
     105    field.project.polyterms[i][0] = (pc11*cd1*field.distort.polyterms[i][0] + pc12*cd2*field.distort.polyterms[i][1]);
     106    field.project.polyterms[i][1] = (pc21*cd1*field.distort.polyterms[i][0] + pc22*cd2*field.distort.polyterms[i][1]);
     107  }
     108  for (i = 0; i < 2; i++) {
     109    field.project.polyterms[0][i] =  field.project.polyterms[0][i] / (cd1*cd1);
     110    field.project.polyterms[1][i] =  field.project.polyterms[1][i] / (cd1*cd2);
     111    field.project.polyterms[2][i] =  field.project.polyterms[2][i] / (cd2*cd2);
     112
     113    field.project.polyterms[3][i] =  field.project.polyterms[3][i] / (cd1*cd1*cd1);
     114    field.project.polyterms[4][i] =  field.project.polyterms[4][i] / (cd1*cd1*cd2);
     115    field.project.polyterms[5][i] =  field.project.polyterms[5][i] / (cd1*cd2*cd2);
     116    field.project.polyterms[6][i] =  field.project.polyterms[6][i] / (cd2*cd2*cd2);
    98117  }
    99118}
    … …  
    116135  }
    117136 
     137  /* separate field into boresite + distortion terms */
     138
    118139  /* bore site center guess */
    119140  strcpy (field.project.ctype, "DEC--TAN");
    … …  
    135156
    136157  /* bore site center guess */
    137   strcpy (field.distort.ctype, "DEC--WRP");
     158  strcpy (field.distort.ctype, "DEC--PLY");
    138159  field.distort.crval1 = 0.0;
    139160  field.distort.crval2 = 0.0;

  • trunk/Ohana/src/mosastro/src/fitpoly.c

    r3323 r3401  
    6868
    6969  int i, j, n, m, M, N;
    70   double max;
     70
     71  if (Npts == 0) {
     72    fprintf (stderr, "warning: no valid pts\n");
     73  }
    7174
    7275  i = 0;
    … …  
    123126
    124127  int i, j, Np, Nv, N;
    125   double c10, c11, c12;
    126   double c20, c21, c22;
    127   double R, D;
     128  double c11, c12;
     129  double c21, c22;
     130  double R;
    128131
    129132  /* update the higher order terms */
    … …  
    163166
    164167  coords[0].Npolyterms = NORDER;
    165   strcpy (coords[0].ctype, "DEC--WRP");
     168  strcpy (coords[0].ctype, "DEC--PLY");
    166169}
    167170
    … …  
    174177
    175178/* NORDER is order of gradient fit : Npolyterms is Norder + 1 */
    176 fit_apply_grads (Coords *distort, Coords *project, int term) {
     179void fit_apply_grads (Coords *distort, Coords *project, int term) {
    177180
    178181  int i, j, Np, Nv1, Nv2;
    179182
    180183  distort[0].Npolyterms = NORDER + 1;
    181   if (distort[0].Npolyterms > 1) strcpy (distort[0].ctype, "DEC--WRP");
     184  if (distort[0].Npolyterms > 1) strcpy (distort[0].ctype, "DEC--PLY");
    182185
    183186  for (i = 0; i < NPOWER + 1; i++) {
    … …  
    212215}
    213216
    214 fit_correct_grads (Gradients *in, Gradients *out, int term) {
     217void fit_correct_grads (Gradients *in, Gradients *out, int term) {
    215218
    216219  int i, k, m, n;

  • trunk/Ohana/src/mosastro/src/gaussj.c

    r3323 r3401  
    1313    ipiv[j] = 0;
    1414
     15  icol = irow = 0;
    1516  for (i = 0; i < n; i++) {
    1617    big = 0.0;

  • trunk/Ohana/src/mosastro/src/greference.c

    r3323 r3401  
    99  if (VERBOSE) fprintf (stderr, "loading astrometric reference data from %s\n", REFCAT);
    1010
     11  stars = NULL;
    1112  Nstars = 0;
    1213  catstats.RA[0]  = field.Rmin;

  • trunk/Ohana/src/mosastro/src/match.c

    r3323 r3401  
    44int match (StarData *refcat, int Nrefcat) {
    55
    6   int i, j, k, N, K, Ntotal;
     6  int i, j, k, K, Ntotal;
    77  int Nmatch, NMATCH;
    88  double dp, dq, radius, Radius;

  • trunk/Ohana/src/mosastro/src/mkmosaic.c

    r3323 r3401  
    33Header *mkmosaic (int Nx, int Ny, int Nstars, Coords *coords) {
    44
    5   int i, Nastro;
     5  int Nastro;
    66  double Xmin, Xmax, Ymin, Ymax, tmp, Cerror, Cprecise;
    77  Header *header;

  • trunk/Ohana/src/mosastro/src/mosastro.c

    r3323 r3401  
    99  StarData *refcat;
    1010  Gradients *grad;
     11
    1112
    1213  ConfigInit (&argc, argv);
    … …  
    6465  if ((DUMP != NULL) && !strcmp (DUMP, "fitchips")) dump_match();
    6566
     67  /* testcoords (); */
    6668  output (argv[2], argv[3]);
    6769  exit (0);

  • trunk/Ohana/src/mosastro/src/output.c

    r3323 r3401  
    2020    strcat (outname, ext);
    2121
     22    /* convert chip from PLY to WRP (linked to field model) */
     23    strcpy (chip[i].map.ctype, "DEC--WRP");
    2224    PutCoords (&chip[i].map, &chip[i].header);
    2325    wchip (outname, &chip[i]);
    … …  
    4143    fclose (f);
    4244  }
     45
     46  {
     47
     48    FILE *f;
     49    double R, D;
     50
     51    fprintf (stderr, "starting test.dat\n");
     52    f = fopen ("test.dat", "w");
     53    if (f == (FILE *) NULL) {
     54      fprintf (stderr, "ERROR: can't create output file \n");
     55      exit (1);
     56    }
     57
     58    RegisterMosaic (&field.project);
     59    for (i = 0; i < chip[0].Nmatch; i++) {
     60      XY_to_RD (&R, &D, chip[0].raw[i].X, chip[0].raw[i].Y, &chip[0].map);
     61      fprintf (f, "%4d  %10.6f %10.6f  %8.2f %8.2f\n",
     62               i, R, D, chip[0].raw[i].X, chip[0].raw[i].Y);
     63    }
     64    fclose (f);
     65  }
     66 
    4367}

  • trunk/Ohana/src/mosastro/src/parse_time.c

    r3294 r3401  
    11# include "mosastro.h"
    22
    3 parse_time (Header *header) {
     3e_time parse_time (Header *header) {
    44
    55  double jd;
    66  int Ny, Nf, mode;
    7   int Nsec, hour, min, sec, year, month, day;
     7  int hour, min, sec, year, month, day;
    88  char *py, *pm, *pd, *c;
    99  char line[256];
     10  e_time Nsec;
    1011
    1112  /* we want to find JD or MJD to get Nsec (seconds since 01/01/1970) */
    … …  
    7576  for (c = strchr (line, 0x2d); c != (char *) NULL; c = strchr (line, 0x2d)) { *c = ' '; }
    7677
     78  Nf = 0;
    7779  switch (mode) {
    7880  case 1:
    … …  
    118120}
    119121
    120 uppercase (char *string) {
     122void uppercase (char *string) {
    121123
    122124  int i;

  • trunk/Ohana/src/mosastro/src/project.c

    r3297 r3401  
    22
    33int project_refcat (StarData *refcat, int Nrefcat) {
    4 
    5   int i;
    64
    75  SkyToTP (refcat, Nrefcat, &field.project);
    … …  
    1311int deproject_stars () {
    1412
    15   int i, j;
     13  int i;
    1614
    1715  for (i = 0; i < Nchip; i++) {
    … …  
    2321int project_stars () {
    2422
    25   int i, j;
     23  int i;
    2624
    2725  for (i = 0; i < Nchip; i++) {
    … …  
    3432int deproject_raw () {
    3533
    36   int i, j;
     34  int i;
    3735
    3836  for (i = 0; i < Nchip; i++) {
    … …  
    4644int project_ref () {
    4745
    48   int i, j;
     46  int i;
    4947
    5048  for (i = 0; i < Nchip; i++) {

  • trunk/Ohana/src/mosastro/src/transforms.c

    r3294 r3401  
    11# include "mosastro.h"
    22
    3 ChipToSky (StarData *stars, int Nstars, Coords *coords) {
     3void ChipToSky (StarData *stars, int Nstars, Coords *coords) {
    44  int i;
    55  for (i = 0; i < Nstars; i++) {
    … …  
    1010}
    1111
    12 ChipToFP (StarData *stars, int Nstars, Coords *coords) {
     12void ChipToFP (StarData *stars, int Nstars, Coords *coords) {
    1313  int i;
    1414  for (i = 0; i < Nstars; i++) {
    … …  
    1717}
    1818
    19 FPtoTP (StarData *stars, int Nstars, Coords *coords) {
     19void FPtoTP (StarData *stars, int Nstars, Coords *coords) {
    2020  int i;
    2121  for (i = 0; i < Nstars; i++) {
    … …  
    2424}
    2525
    26 TPtoSky (StarData *stars, int Nstars, Coords *coords) {
     26void TPtoSky (StarData *stars, int Nstars, Coords *coords) {
    2727  int i;
    2828  for (i = 0; i < Nstars; i++) {
    … …  
    3333}
    3434
    35 SkyToTP (StarData *stars, int Nstars, Coords *coords) {
     35void SkyToTP (StarData *stars, int Nstars, Coords *coords) {
    3636  int i;
    3737  for (i = 0; i < Nstars; i++) {
    … …  
    4040}
    4141
    42 TPtoFP (StarData *stars, int Nstars, Coords *coords) {
     42void TPtoFP (StarData *stars, int Nstars, Coords *coords) {
    4343  int i;
    4444  for (i = 0; i < Nstars; i++) {
    … …  
    4747}
    4848
    49 FPtoChip (StarData *stars, int Nstars, Coords *coords) {
     49void FPtoChip (StarData *stars, int Nstars, Coords *coords) {
    5050  int i;
    5151  for (i = 0; i < Nstars; i++) {
Note: See TracChangeset for help on using the changeset viewer.