IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Ignore:
Timestamp:
Dec 3, 2015, 4:16:29 PM (11 years ago)
Author:
eugene
Message:

merge changes from eam_branches/ipp-20151113

Location:
trunk
Files:
10 edited

Legend:

Unmodified
Added
Removed
  • trunk

  • trunk/Ohana

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

    r38986 r39225  
    254254int SRC_MEAS_TOOFEW; //catalog objects wich fewer detections then this are ignored
    255255double MIN_ERROR;
     256
     257double MIN_DISTANCE_MOD; 
     258double MAX_DISTANCE_MOD; 
     259double MAX_DISTANCE_MOD_ERR;
    256260
    257261int    IMFIT_TOO_FEW; // need more than this number of stars to fit an image
     
    323327PhotCode     *DCR_RED_PHOTCODE_POS, *DCR_RED_PHOTCODE_NEG;
    324328int           DCR_RED_NSEC_POS,      DCR_RED_NSEC_NEG;
     329
     330float *LoopWeight2MASS;
     331float *LoopWeightTycho;
    325332
    326333int ImagSelect;
  • trunk/Ohana/src/relastro/src/BrightCatalog.c

    r38553 r39225  
    159159    GET_COLUMN(uR,             "U_RA",           float);
    160160    GET_COLUMN(uD,             "U_DEC",          float);
    161     GET_COLUMN(duR,            "V_RA_ERR",       float);
    162     GET_COLUMN(duD,            "V_DEC_ERR",      float);
     161    GET_COLUMN(duR,            "U_RA_ERR",       float);
     162    GET_COLUMN(duD,            "U_DEC_ERR",      float);
     163    GET_COLUMN(uRgal,          "U_RA_GAL",       float);
     164    GET_COLUMN(uDgal,          "U_DEC_GAL",      float);
    163165    GET_COLUMN(P,              "PAR",            float);
    164166    GET_COLUMN(dP,             "PAR_ERR",        float);
     
    189191    ALLOCATE (average, Average, Nrow);
    190192    for (i = 0; i < Nrow; i++) {
     193      dvo_average_init (&average[i]);
    191194      average[i].R               = R[i]               ;
    192195      average[i].D               = D[i]               ;
     
    197200      average[i].duR             = duR[i]             ;
    198201      average[i].duD             = duD[i]             ;
     202      average[i].uRgal           = uRgal[i]           ;
     203      average[i].uDgal           = uDgal[i]           ;
    199204      average[i].P               = P[i]               ;
    200205      average[i].dP              = dP[i]              ;
     
    230235    free (duR);
    231236    free (duD);
     237    free (uRgal);
     238    free (uDgal);
    232239    free (P);
    233240    free (dP);
     
    489496    gfits_define_bintable_column (&theader, "E", "U_RA",           "RA*cos(D) proper-motion",                             "", 1.0, 0.0);
    490497    gfits_define_bintable_column (&theader, "E", "U_DEC",          "DEC proper-motion",                                   "", 1.0, 0.0);
    491     gfits_define_bintable_column (&theader, "E", "V_RA_ERR",       "RA*cos(D) p-m error",                                 "", 1.0, 0.0);
    492     gfits_define_bintable_column (&theader, "E", "V_DEC_ERR",      "DEC p-m error",                                       "", 1.0, 0.0);
     498    gfits_define_bintable_column (&theader, "E", "U_RA_ERR",       "RA*cos(D) p-m error",                                 "", 1.0, 0.0);
     499    gfits_define_bintable_column (&theader, "E", "U_DEC_ERR",      "DEC p-m error",                                       "", 1.0, 0.0);
     500    gfits_define_bintable_column (&theader, "E", "U_RA_GAL",       "RA*cos(D) p-m error",                                 "", 1.0, 0.0);
     501    gfits_define_bintable_column (&theader, "E", "U_DEC_GAL",      "DEC p-m error",                                       "", 1.0, 0.0);
    493502    gfits_define_bintable_column (&theader, "E", "PAR",            "parallax",                                            "", 1.0, 0.0);
    494503    gfits_define_bintable_column (&theader, "E", "PAR_ERR",        "parallax error",                                      "", 1.0, 0.0);
     
    526535    float    *duR           ; ALLOCATE (duR           , float   , catalog->Naverage);
    527536    float    *duD           ; ALLOCATE (duD           , float   , catalog->Naverage);
     537    float    *uRgal         ; ALLOCATE (uRgal         , float   , catalog->Naverage);
     538    float    *uDgal         ; ALLOCATE (uDgal         , float   , catalog->Naverage);
    528539    float    *P             ; ALLOCATE (P             , float   , catalog->Naverage);
    529540    float    *dP            ; ALLOCATE (dP            , float   , catalog->Naverage);
     
    560571      duR[i]             = average[i].duR             ;
    561572      duD[i]             = average[i].duD             ;
     573      uRgal[i]           = average[i].uRgal           ;
     574      uDgal[i]           = average[i].uDgal           ;
    562575      P[i]               = average[i].P               ;
    563576      dP[i]              = average[i].dP              ;
     
    591604    gfits_set_bintable_column (&theader, &ftable, "U_RA",           uR,              catalog->Naverage);
    592605    gfits_set_bintable_column (&theader, &ftable, "U_DEC",          uD,              catalog->Naverage);
    593     gfits_set_bintable_column (&theader, &ftable, "V_RA_ERR",       duR,             catalog->Naverage);
    594     gfits_set_bintable_column (&theader, &ftable, "V_DEC_ERR",      duD,             catalog->Naverage);
     606    gfits_set_bintable_column (&theader, &ftable, "U_RA_ERR",       duR,             catalog->Naverage);
     607    gfits_set_bintable_column (&theader, &ftable, "U_DEC_ERR",      duD,             catalog->Naverage);
     608    gfits_set_bintable_column (&theader, &ftable, "U_RA_GAL",       uRgal,           catalog->Naverage);
     609    gfits_set_bintable_column (&theader, &ftable, "U_DEC_GAL",      uDgal,           catalog->Naverage);
    595610    gfits_set_bintable_column (&theader, &ftable, "PAR",            P,               catalog->Naverage);
    596611    gfits_set_bintable_column (&theader, &ftable, "PAR_ERR",        dP,              catalog->Naverage);
     
    604619    gfits_set_bintable_column (&theader, &ftable, "NMEASURE",       Nmeasure,        catalog->Naverage);
    605620    gfits_set_bintable_column (&theader, &ftable, "NMISSING",       Nmissing,        catalog->Naverage);
    606     gfits_set_bintable_column (&theader, &ftable, "NGALPHOT",      Ngalphot,       catalog->Naverage);
     621    gfits_set_bintable_column (&theader, &ftable, "NGALPHOT",       Ngalphot,        catalog->Naverage);
    607622    gfits_set_bintable_column (&theader, &ftable, "OFF_MEASURE",    measureOffset,   catalog->Naverage);
    608623    gfits_set_bintable_column (&theader, &ftable, "OFF_MISSING",    missingOffset,   catalog->Naverage);
     
    624639    free (duR);
    625640    free (duD);
     641    free (uRgal);
     642    free (uDgal);
    626643    free (P);
    627644    free (dP);
  • trunk/Ohana/src/relastro/src/ConfigInit.c

    r37807 r39225  
    3535  if (!ScanConfig (config, "RELASTRO_DPOS_MAX",      "%lf", 0, &DPOS_MAX))         DPOS_MAX = 6.0;   
    3636  if (!ScanConfig (config, "ADDSTAR_RADIUS",         "%lf", 0, &ADDSTAR_RADIUS))   ADDSTAR_RADIUS = 1.0;
     37
     38  if (!ScanConfig (config, "RELASTRO_MIN_DISTANCE_MOD",     "%lf", 0, &MIN_DISTANCE_MOD))     MIN_DISTANCE_MOD     =  7.5;
     39  if (!ScanConfig (config, "RELASTRO_MAX_DISTANCE_MOD",     "%lf", 0, &MAX_DISTANCE_MOD))     MAX_DISTANCE_MOD     = 15.0;
     40  if (!ScanConfig (config, "RELASTRO_MAX_DISTANCE_MOD_ERR", "%lf", 0, &MAX_DISTANCE_MOD_ERR)) MAX_DISTANCE_MOD_ERR =  1.0;
    3741
    3842  if (!ScanConfig (config, "USE_FIXED_PIXCOORDS", "%d", 0, &USE_FIXED_PIXCOORDS))  USE_FIXED_PIXCOORDS = FALSE;
  • trunk/Ohana/src/relastro/src/GetAstromError.c

    r38986 r39225  
    7373  int has2MASS = USE_GALAXY_MODEL &&  isImage && (measure[0].dbFlags & ID_MEAS_OBJECT_HAS_2MASS);
    7474
    75   // serious hack: if the object has 2MASS, we set this internal bit and adjust the
     75  // modest hack: if the object has 2MASS or Tycho, we set this internal bit and adjust the
    7676  // weight to ensure the image is tied down to the 2mass frame
    7777
    78   if (has2MASS) {
    79     switch (Nloop) {
    80       case 0:
    81         dPtotal = dPtotal / 1000.0;
    82         break;
    83       case 1:
    84         dPtotal = dPtotal / 300.0;
    85         break;
    86       case 2:
    87         dPtotal = dPtotal / 100.0;
    88         break;
    89       case 3:
    90         dPtotal = dPtotal / 30.0;
    91         break;
    92       case 4:
    93         dPtotal = dPtotal / 10.0;
    94         break;
    95       case 5:
    96         dPtotal = dPtotal / 10.0;
    97         break;
    98       default:
    99         break;
    100     }
     78  if (has2MASS && LoopWeight2MASS) {
     79    dPtotal = dPtotal / LoopWeight2MASS[Nloop];
    10180  }
    102   if (is2MASS) {
    103     switch (Nloop) {
    104       case 0:
    105         dPtotal = dPtotal / 1000.0;
    106         break;
    107       case 1:
    108         dPtotal = dPtotal / 300.0;
    109         break;
    110       case 2:
    111         dPtotal = dPtotal / 100.0;
    112         break;
    113       case 3:
    114         dPtotal = dPtotal / 30.0;
    115         break;
    116       case 4:
    117         dPtotal = dPtotal / 10.0;
    118         break;
    119       case 5:
    120         dPtotal = dPtotal / 10.0;
    121         break;
    122       default:
    123         break;
    124     }
     81  if (is2MASS && LoopWeight2MASS) {
     82    dPtotal = dPtotal / LoopWeight2MASS[Nloop];
    12583  }
    126   if (isTycho) {
    127     switch (Nloop) {
    128       case 0:
    129       case 1:
    130         dPtotal = dPtotal / 200.0;
    131         break;
    132       case 2:
    133       case 3:
    134         dPtotal = dPtotal / 100.0;
    135         break;
    136       case 4:
    137       case 5:
    138         dPtotal = dPtotal / 50.0;
    139         break;
    140       default:
    141         break;
    142     }
     84  if (isTycho && LoopWeightTycho) {
     85    dPtotal = dPtotal / LoopWeightTycho[Nloop];
    14386  }
    144 
    14587  return (dPtotal);
    14688}
     
    197139
    198140  dPtotal = MAX (dPtotal, MIN_ERROR);
     141
     142  // early on, we want 2MASS and Tycho to have a very high weight.  This will force images
     143  // to match the 2MASS / Tycho / ICRS reference frame.  As Nloop gets higher, the weight
     144  // needs to drop to allow the ps1 measurements to drive the solution
     145  int is2MASS  = USE_GALAXY_MODEL && !isImage && (measure[0].photcode >= 2011) && (measure[0].photcode <= 2013);
     146  int isTycho  = USE_GALAXY_MODEL && !isImage && (measure[0].photcode >= 2020) && (measure[0].photcode <= 2021);
     147  int has2MASS = USE_GALAXY_MODEL &&  isImage && (measure[0].dbFlags & ID_MEAS_OBJECT_HAS_2MASS);
     148
     149  // modest hack: if the object has 2MASS or Tycho, we set this internal bit and adjust the
     150  // weight to ensure the image is tied down to the 2mass frame
     151
     152  if (has2MASS && LoopWeight2MASS) {
     153    dPtotal = dPtotal / LoopWeight2MASS[Nloop];
     154  }
     155  if (is2MASS && LoopWeight2MASS) {
     156    dPtotal = dPtotal / LoopWeight2MASS[Nloop];
     157  }
     158  if (isTycho && LoopWeightTycho) {
     159    dPtotal = dPtotal / LoopWeightTycho[Nloop];
     160  }
     161
    199162  return (dPtotal);
    200163}
    201164
    202 /* for a long time, psphot was either not reported position errors, or was reporting
     165/* for a long time, psphot was either not reporting position errors, or was reporting
    203166 * completely wrong astrometry errors.  This function lets us handle, in the
    204167 * configuration, different strategies to generating a position error
  • trunk/Ohana/src/relastro/src/args.c

    r38986 r39225  
    44void usage_merge_source (void);
    55void usage_merge_source_id (char *name);
     6float *ParseLoopWeights (char *rawlist);
    67
    78int args (int argc, char **argv) {
     
    529530    remove_argument (N, &argc, argv);
    530531    NLOOP = atof (argv[N]);
     532    remove_argument (N, &argc, argv);
     533  }
     534
     535  // e.g., -loop-weights-2mass 1000,300,300,200,200,100
     536  // NOTE: this must come after -nloop above
     537  LoopWeight2MASS = NULL;
     538  if ((N = get_argument (argc, argv, "-loop-weights-2mass"))) {
     539    remove_argument (N, &argc, argv);
     540    LoopWeight2MASS = ParseLoopWeights (argv[N]);
     541    remove_argument (N, &argc, argv);
     542  }
     543  LoopWeightTycho = NULL;
     544  if ((N = get_argument (argc, argv, "-loop-weights-tycho"))) {
     545    remove_argument (N, &argc, argv);
     546    LoopWeightTycho = ParseLoopWeights (argv[N]);
    531547    remove_argument (N, &argc, argv);
    532548  }
     
    9991015}
    10001016
     1017float *ParseLoopWeights (char *rawlist) {
     1018
     1019  float *weights = NULL;
     1020  ALLOCATE (weights, float, NLOOP);
     1021
     1022  int Nloop = 0;
     1023
     1024  /* parse the comma-separated list of photcodes */
     1025  char *myList = strcreate(rawlist);
     1026  char *list = myList;
     1027  char *entry = NULL;
     1028  char *ptr = NULL;
     1029  while ((Nloop < NLOOP) && ((entry = strtok_r (list, ",", &ptr)) != NULL)) {
     1030    list = NULL; // pass NULL on successive strtok_r calls
     1031
     1032    weights[Nloop] = atof(entry);
     1033    Nloop ++;
     1034  }
     1035  free (myList);
     1036
     1037  if (Nloop == 0) {
     1038    fprintf (stderr, "syntax error parsing weights: %s\n", rawlist);
     1039    exit (3);
     1040  }
     1041
     1042  while (Nloop < NLOOP) {
     1043    weights[Nloop] = weights[Nloop - 1];
     1044    Nloop ++;
     1045  }
     1046  return weights;
     1047}
  • trunk/Ohana/src/relastro/src/bcatalog.c

    r38986 r39225  
    117117      if (!isfinite(catalog[0].average[i].uRgal)) continue;
    118118      if (!isfinite(catalog[0].average[i].uDgal)) continue;
     119     
     120      // filter stars which are not well-handled by the galactic motion model
     121      int mStarpar = catalog[0].average[i].starparOffset;
     122      myAssert (mStarpar < catalog[0].Nstarpar, "oops");
     123     
     124      StarPar *starpar = &catalog[0].starpar[mStarpar];
     125      if (starpar->DistMag  < MIN_DISTANCE_MOD) continue;
     126      if (starpar->DistMag  > MAX_DISTANCE_MOD) continue;
     127      if (starpar->dDistMag > MAX_DISTANCE_MOD_ERR) continue;
    119128    }
    120129
  • trunk/Ohana/src/relastro/src/launch_region_hosts.c

    r38986 r39225  
    134134    if (MaxMeanOffset != 10.0) strextend (&command, "-max-mean-offset  %f", MaxMeanOffset);
    135135
     136    if (LoopWeight2MASS) {
     137      char *weightline = NULL;
     138      for (i = 0; i < NLOOP; i++) {
     139        strextend (&weightline, "%f,", LoopWeight2MASS[i]);
     140      }
     141      strextend (&command, "-loop-weights-2mass %s", weightline);
     142      free (weightline);
     143    }
     144    if (LoopWeightTycho) {
     145      char *weightline = NULL;
     146      for (i = 0; i < NLOOP; i++) {
     147        strextend (&weightline, "%f,", LoopWeightTycho[i]);
     148      }
     149      strextend (&command, "-loop-weights-tycho %s", weightline);
     150      free (weightline);
     151    }
     152
    136153    strextend (&command, "-D RELASTRO_SRC_MEAS_TOOFEW %d", SRC_MEAS_TOOFEW);
    137154    strextend (&command, "-D RELASTRO_SIGMA_LIM %f", SIGMA_LIM);
  • trunk/Ohana/src/relastro/src/load_catalogs.c

    r38986 r39225  
    4444    pcatalog->catflags  = DVO_LOAD_AVERAGE | DVO_LOAD_MEASURE | DVO_LOAD_SECFILT; // don't need to load all data at this point
    4545    pcatalog->Nsecfilt  = GetPhotcodeNsecfilt ();
     46
     47    if (USE_GALAXY_MODEL) {
     48      pcatalog->catflags = pcatalog->catflags | DVO_LOAD_STARPAR;
     49    }
    4650
    4751    // loads Average, Measure, SecFilt
     
    171175    strextend (&command, " -minerror %f", MIN_ERROR);
    172176    strextend (&command, " -D RELASTRO_SIGMA_LIM %f", SIGMA_LIM);
     177    strextend (&command, " -D RELASTRO_MIN_DISTANCE_MOD %f",     MIN_DISTANCE_MOD);
     178    strextend (&command, " -D RELASTRO_MAX_DISTANCE_MOD %f",     MAX_DISTANCE_MOD);
     179    strextend (&command, " -D RELASTRO_MAX_DISTANCE_MOD_ERR %f", MAX_DISTANCE_MOD_ERR);
    173180
    174181    if (FIT_MODE == FIT_PM_ONLY)         strextend (&command, "-pm");
Note: See TracChangeset for help on using the changeset viewer.