IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 37606


Ignore:
Timestamp:
Nov 16, 2014, 1:35:31 PM (12 years ago)
Author:
eugene
Message:

galaxy model correction basically works (need to be careful about starting position)

Location:
branches/eam_branches/ipp-20140904/Ohana/src/relastro
Files:
13 edited

Legend:

Unmodified
Added
Removed

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/include/relastro.h

    r37604 r37606  
    181181double MIN_ERROR;
    182182
     183int    IMFIT_TOO_FEW; // need more than this number of stars to fit an image
    183184int    IMFIT_CLIP_NITER; // number of clipping iterations to perform in FitChip
    184185double IMFIT_CLIP_NSIGMA; // number of sigma to clip in FitChip
    … …  
    409410void FitSimple (StarData *raw, StarData *ref, int Nmatch, Coords *coords);
    410411void initObjectData (Catalog *catalog, int Ncatalog);
    411 int UpdateObjects (Catalog *catalog, int Ncatalog);
     412int UpdateObjects (Catalog *catalog, int Ncatalog, int Nloop);
    412413int UpdateSimple (Catalog *catalog, int Ncatalog);
    413414int UpdateChips (Catalog *catalog, int Ncatalog);

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ConfigInit.c

    r37545 r37606  
    2525  if (!ScanConfig (config, "RELASTRO_IMFIT_CLIP_NSIGMA",   "%lf", 0, &IMFIT_CLIP_NSIGMA))   IMFIT_CLIP_NSIGMA = 3.0;
    2626  if (!ScanConfig (config, "RELASTRO_IMFIT_SYS_SIGMA_LIM", "%lf", 0, &IMFIT_SYS_SIGMA_LIM)) IMFIT_SYS_SIGMA_LIM = 0.01;
     27  if (!ScanConfig (config, "RELASTRO_IMFIT_TOO_FEW",       "%d",  0, &IMFIT_TOO_FEW))       IMFIT_TOO_FEW = 5;
    2728
    2829  if (!ScanConfig (config, "PM_DT_MIN",              "%lf", 0, &PM_DT_MIN))        PM_DT_MIN = 0.25;   

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/FitChip.c

    r37604 r37606  
    5555    // residual vector field
    5656    // dx vs mag, dy vs mag
     57    if (relastroGetVisual()) {
     58      fprintf (stderr, "fitting chip %s, %d matched stars\n", image->name, Nmatch);
     59    }
    5760    relastroVisualPlotChipFit(raw, ref, dRmax, Nmatch);
    5861
    … …  
    7477    // XXX WATCH OUT HERE!!
    7578    if (!CHIPMAP) order_use = MIN(order_use, 1);
    76     fprintf (stderr, "using %d for %s\n", order_use, image[0].name);
     79    if (VERBOSE2) fprintf (stderr, "using %d for %s\n", order_use, image[0].name);
    7780
    7881    // when fitting the map, first fit a linear model (below? change Npolyterms to -1)

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/ImageOps.c

    r37604 r37606  
    491491
    492492    int TESTPT = FALSE;
    493 
    494     TESTPT |= CAT_ID_SRC && OBJ_ID_SRC && (measure->imageID == CAT_ID_SRC) && (measure->detID == OBJ_ID_SRC);
    495     TESTPT |= CAT_ID_DST && OBJ_ID_DST && (measure->imageID == CAT_ID_DST) && (measure->detID == OBJ_ID_DST);
     493    TESTPT |= CAT_ID_SRC && OBJ_ID_SRC && (measure->catID == CAT_ID_SRC) && (measure->objID == OBJ_ID_SRC);
     494    TESTPT |= CAT_ID_DST && OBJ_ID_DST && (measure->catID == CAT_ID_DST) && (measure->objID == OBJ_ID_DST);
    496495    if (TESTPT) {
    497496      fprintf (stderr, "got test det\n");
    … …  
    736735
    737736    n = measure[0].averef;
     737
     738    int TESTPT = FALSE;
     739    TESTPT |= CAT_ID_SRC && OBJ_ID_SRC && (catalog[c].average[n].catID == CAT_ID_SRC) && (catalog[c].average[n].objID == OBJ_ID_SRC);
     740    TESTPT |= CAT_ID_DST && OBJ_ID_DST && (catalog[c].average[n].catID == CAT_ID_DST) && (catalog[c].average[n].objID == OBJ_ID_DST);
     741    if (TESTPT) {
     742      fprintf (stderr, "got test det\n");
     743    }
    738744
    739745    // an object with only one detection provides no information about the image calibration

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/UpdateChips.c

    r37562 r37606  
    7474      Nskip ++;
    7575      mode[i] = 0;
     76      continue;
     77    }
     78    if (Nraw <= IMFIT_TOO_FEW) {
     79      Nskip ++;
     80      mode[i] = 0;
     81      free (raw);
    7682      continue;
    7783    }

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/UpdateObjectOffsets.c

    r37038 r37606  
    5959
    6060    // uses MeasureTiny values, but will also update Measure values when populated
    61     UpdateObjects (&catalog, 1);
     61    UpdateObjects (&catalog, 1, 0);
    6262
    6363    free_tiny_values(&catalog);

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/UpdateObjects.c

    r37604 r37606  
    22# define PAR_TOOFEW 5
    33
    4 int UpdateObjects_Chips (Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int Nsecfilt, FitStats *fitStats, int i, off_t m);
     4int UpdateObjects_Chips (Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int Nsecfilt, FitStats *fitStats, int i, off_t m, int Nloop);
    55int UpdateObjects_Stack (Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int Nsecfilt, FitStats *fitStats);
    66
    … …  
    9999// This function operates on both Measure and MeasureTiny.  In the big stages, this should
    100100// be called with just MeasureTiny set and Measure == NULL
    101 int UpdateObjects (Catalog *catalog, int Ncatalog) {
     101int UpdateObjects (Catalog *catalog, int Ncatalog, int Nloop) {
    102102
    103103  initObjectData (catalog, Ncatalog);
    … …  
    130130
    131131      UpdateObjects_Stack(average, secfilt, measure, measureBig, Nsecfilt, &fitStatsStack);
    132       UpdateObjects_Chips(average, secfilt, measure, measureBig, Nsecfilt, &fitStatsChips, i, m);
     132      UpdateObjects_Chips(average, secfilt, measure, measureBig, Nsecfilt, &fitStatsChips, i, m, Nloop);
    133133    }
    134134    if (VERBOSE) fprintf (stderr, "catalog %d : chips "OFF_T_FMT" ave, "OFF_T_FMT" pm, "OFF_T_FMT" par : Nskip "OFF_T_FMT", Noffset "OFF_T_FMT"\n",  i,  fitStatsChips.Nave,  fitStatsChips.Npm,  fitStatsChips.Npar,  fitStatsChips.Nskip, fitStatsChips.Noffset);
    … …  
    146146// This function operates on both Measure and MeasureTiny.  In the big stages, this should
    147147// be called with just MeasureTiny set and Measure == NULL
    148 int UpdateObjects_Chips (Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int Nsecfilt, FitStats *fitStats, int i, off_t m) {
     148int UpdateObjects_Chips (Average *average, SecFilt *secfilt, MeasureTiny *measure, Measure *measureBig, int Nsecfilt, FitStats *fitStats, int i, off_t m, int Nloop) {
    149149
    150150  int setRefColor = areImagesMatched();
    … …  
    209209    double Di = getMeanD (&measure[k], average, secfilt);
    210210
    211     if (USE_GALAXY_MODEL) {
     211    if (USE_GALAXY_MODEL && (Nloop > 0)) {
    212212      Ri -= measure[k].RoffGAL / 3600.0;
    213213      Di -= measure[k].DoffGAL / 3600.0;
    … …  
    470470      break;
    471471  }
    472   if (XVERB) fprintf (stderr, "%f %f -> %f %f (%f,%f) pm=(%f %f) plx=(%f +/- %f)\n",
    473                       average[0].R,
    474                       average[0].D,
    475                       fit.Ro, fit.Do,
    476                       3600*(average[0].R - fit.Ro),
    477                       3600*(average[0].D - fit.Do),
    478                       fit.uR, fit.uD, fit.p, fit.dp);
     472  if (XVERB) {
     473    fprintf (stderr, "%f %f -> %f %f (%f,%f) pm=(%f %f) plx=(%f +/- %f)\n",
     474             average[0].R,
     475             average[0].D,
     476             fit.Ro, fit.Do,
     477             3600*(average[0].R - fit.Ro),
     478             3600*(average[0].D - fit.Do),
     479             fit.uR, fit.uD, fit.p, fit.dp);
     480    fprintf (stderr, "-----\n");
     481  }
    479482
    480483  // make sure that the fit succeeded
    … …  
    510513  }
    511514
    512   average[0].R          = fit.Ro; // RA in degrees
    513   average[0].D          = fit.Do; // DEC in degrees
    514   average[0].dR         = fit.dRo; // RA scatter in arcsec
    515   average[0].dD         = fit.dDo; // DEC scatter in arcsec
    516 
    517   average[0].uR         = fit.uR; // RA proper motion in arcsec/year
    518   average[0].uD         = fit.uD; // DEC proper motion in arcsec/year
    519   average[0].duR        = fit.duR; // RA proper motion error in arcsec/year
    520   average[0].duD        = fit.duD; // DEC proper motion error in arcsec/year
    521 
    522   average[0].P          = fit.p; // parallax in arcsec
    523   average[0].dP         = fit.dp; // parallax error in arcsec
    524 
    525   average[0].ChiSqAve   = fitAve.chisq;
    526   average[0].ChiSqPM    = fitPM.chisq;
    527   average[0].ChiSqPar   = fitPAR.chisq;
    528 
    529   average[0].Tmean      = (Tmean * 86400 * 365.25) + T2000;
    530   average[0].Trange     = (Trange * 86400 * 365.25);
    531   average[0].Npos       = fit.Nfit;
    532 
    533   // unset the NO_ASTROM bit (not(NO_ASTROM) == HAVE_ASTROM)
    534   average[0].flags &= ~ID_STAR_NO_ASTROM;
    535 
    536515  if (XVERB) fprintf (stderr, "%f %f -> %f %f (%f,%f) pm=(%f %f) chisq=(%f, %f, %f)\n",
    537516                      average[0].R,
    … …  
    544523                      fitAve.chisq, fitPM.chisq, fitPAR.chisq);
    545524
     525  average[0].R          = fit.Ro; // RA in degrees
     526  average[0].D          = fit.Do; // DEC in degrees
     527  average[0].dR         = fit.dRo; // RA scatter in arcsec
     528  average[0].dD         = fit.dDo; // DEC scatter in arcsec
     529
     530  average[0].uR         = fit.uR; // RA proper motion in arcsec/year
     531  average[0].uD         = fit.uD; // DEC proper motion in arcsec/year
     532  average[0].duR        = fit.duR; // RA proper motion error in arcsec/year
     533  average[0].duD        = fit.duD; // DEC proper motion error in arcsec/year
     534
     535  average[0].P          = fit.p; // parallax in arcsec
     536  average[0].dP         = fit.dp; // parallax error in arcsec
     537
     538  average[0].ChiSqAve   = fitAve.chisq;
     539  average[0].ChiSqPM    = fitPM.chisq;
     540  average[0].ChiSqPar   = fitPAR.chisq;
     541
     542  average[0].Tmean      = (Tmean * 86400 * 365.25) + T2000;
     543  average[0].Trange     = (Trange * 86400 * 365.25);
     544  average[0].Npos       = fit.Nfit;
     545
     546  // unset the NO_ASTROM bit (not(NO_ASTROM) == HAVE_ASTROM)
     547  average[0].flags &= ~ID_STAR_NO_ASTROM;
     548
    546549  return (TRUE);
    547550}

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/high_speed_objects.c

    r37386 r37606  
    266266
    267267  populate_tiny_values (&catalogOut, DVO_TV_MEASURE);
    268   UpdateObjects (&catalogOut, 1);
     268  UpdateObjects (&catalogOut, 1, 0);
    269269  free_tiny_values(&catalogOut);
    270270

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/hpm_objects.c

    r37386 r37606  
    252252    // distinguishing a bad fit from a fit with 1 or 2 bad points
    253253    FIT_MODE = FIT_PM_ONLY;
    254     UpdateObjects (&testcat, 1);
     254    UpdateObjects (&testcat, 1, 0);
    255255
    256256    free_tiny_values(&testcat);
    … …  
    282282
    283283  populate_tiny_values (&catalogOut, DVO_TV_MEASURE);
    284   UpdateObjects (&catalogOut, 1);
     284  UpdateObjects (&catalogOut, 1, 0);
    285285  free_tiny_values(&catalogOut);
    286286

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/relastroVisual.c

    r37038 r37606  
    2424/** Spawn a kapa window */
    2525static int initWindow(int *kapid) {
     26  if (NTHREADS) {
     27    fprintf (stderr, "-visual and -threads are incompatible, quitting\n");
     28    exit (2);
     29  }
    2630  if (*kapid == -1) {
    2731    *kapid = KapaOpenNamedSocket("kapa", "relastro_plots");

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/relastro_images.c

    r37545 r37606  
    6565    case TARGET_SIMPLE:
    6666      for (i = 0; i < NLOOP; i++) {
    67         UpdateObjects (catalog, Ncatalog);
     67        UpdateObjects (catalog, Ncatalog, i);
    6868        UpdateSimple (catalog, Ncatalog);
    6969      }
    … …  
    7373      for (i = 0; i < NLOOP; i++) {
    7474        // XXX is is really possible that I do not update the measure.R,D on each loop??
    75         UpdateObjects (catalog, Ncatalog); // calculate <R>,<D>
    76         FrameCorrection (catalog, Ncatalog);
     75        UpdateObjects (catalog, Ncatalog, i); // calculate <R>,<D>
     76        if (i > 2) FrameCorrection (catalog, Ncatalog);
    7777        UpdateChips (catalog, Ncatalog);   // measure.X,Y -> R,D, fit image.coords
    7878        MARKTIME("update chips: %f sec\n", dtime);
    … …  
    8484    case TARGET_MOSAICS:
    8585      for (i = 0; i < NLOOP; i++) {
    86         UpdateObjects (catalog, Ncatalog);
     86        UpdateObjects (catalog, Ncatalog, i);
    8787        UpdateMosaic (catalog, Ncatalog);
    8888      }

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/relastro_objects.c

    r36833 r37606  
    7777    populate_tiny_values(&catalog, DVO_TV_MEASURE);
    7878
    79     UpdateObjects (&catalog, 1);
     79    UpdateObjects (&catalog, 1, 0);
    8080
    8181    free_tiny_values(&catalog);

  • branches/eam_branches/ipp-20140904/Ohana/src/relastro/src/relastro_parallel_images.c

    r36833 r37606  
    7676    case TARGET_SIMPLE:
    7777      for (i = 0; i < NLOOP; i++) {
    78         UpdateObjects (catalog, Ncatalog);
     78        UpdateObjects (catalog, Ncatalog, i);
    7979        share_mean_pos (catalog, Ncatalog, regionHosts, i);
    8080        slurp_mean_pos (catalog, Ncatalog, regionHosts, i);
    … …  
    8787    case TARGET_CHIPS:
    8888      for (i = 0; i < NLOOP; i++) {
    89         UpdateObjects (catalog, Ncatalog);
     89        UpdateObjects (catalog, Ncatalog, i);
    9090        share_mean_pos (catalog, Ncatalog, regionHosts, i);
    9191        slurp_mean_pos (catalog, Ncatalog, regionHosts, i);
    … …  
    101101    case TARGET_MOSAICS:
    102102      for (i = 0; i < NLOOP; i++) {
    103         UpdateObjects (catalog, Ncatalog);
     103        UpdateObjects (catalog, Ncatalog, i);
    104104        share_mean_pos (catalog, Ncatalog, regionHosts, i);
    105105        slurp_mean_pos (catalog, Ncatalog, regionHosts, i);
Note: See TracChangeset for help on using the changeset viewer.