Index: trunk/psModules/src/objects/pmPSFtryMetric.c =================================================================== --- trunk/psModules/src/objects/pmPSFtryMetric.c (revision 26190) +++ trunk/psModules/src/objects/pmPSFtryMetric.c (revision 26260) @@ -76,4 +76,5 @@ pmSourceVisualPlotPSFMetric (psfTry); + pmSourceVisualPlotPSFMetricSubpix (psfTry); psFree (stats); Index: trunk/psModules/src/objects/pmPSFtryModel.c =================================================================== --- trunk/psModules/src/objects/pmPSFtryModel.c (revision 26190) +++ trunk/psModules/src/objects/pmPSFtryModel.c (revision 26260) @@ -72,9 +72,5 @@ } - // XXX set the min number of needed source more carefully (depends on psfTrendMode?) - int orderMax = PS_MAX (options->psfTrendNx, options->psfTrendNy); - if ((sources->n < 15) && (orderMax >= 3)) orderMax = 2; - if ((sources->n < 11) && (orderMax >= 2)) orderMax = 1; - if ((sources->n < 8) && (orderMax >= 1)) orderMax = 0; + // hard limit on minimum number of stars if ((sources->n < 3)) { psError (PS_ERR_UNKNOWN, true, "failed to determine PSF parameters"); @@ -82,11 +78,23 @@ } - int orderMin; + // this is a bit tricky, because we have two cases (MAP vs POLY), and they have a different + // definition for 'order' (order_MAP = order_POLY + 1). in addition, we have a + // user-specified MAX order, which we should respect, regardless of the mode + + // set the max order (0 = constant) which the number of psf stars can support: + // rule of thumb: require 3 stars per 'cell' (order+1)^2 + int MaxOrderForStars = 0; + if (sources->n >= 12) MaxOrderForStars = 1; // 4 cells + if (sources->n >= 27) MaxOrderForStars = 2; // 9 cells + if (sources->n >= 48) MaxOrderForStars = 3; // 16 cells + if (sources->n > 75) MaxOrderForStars = 4; // 25 cells + + int orderMax = PS_MAX (options->psfTrendNx, options->psfTrendNy); + int orderMin = 0; if (options->psfTrendMode == PM_TREND_MAP) { - orderMin = 1; - orderMax = PS_MAX(orderMax, 1); - } else { - orderMin = 0; - } + MaxOrderForStars ++; + orderMin ++; + } + orderMax = PS_MIN (orderMax, MaxOrderForStars); // save the raw source mask (generated by pmPSFtryFitEXT) Index: trunk/psModules/src/objects/pmSourceIO.c =================================================================== --- trunk/psModules/src/objects/pmSourceIO.c (revision 26190) +++ trunk/psModules/src/objects/pmSourceIO.c (revision 26260) @@ -980,4 +980,11 @@ if (!tableHeader) psAbort("cannot read table header"); + char *xtension = psMetadataLookupStr (NULL, tableHeader, "XTENSION"); + if (!xtension) psAbort("cannot read table type"); + if (strcmp (xtension, "BINTABLE")) { + psWarning ("no binary table in extension %s, skipping\n", dataname); + return false; + } + char *exttype = psMetadataLookupStr (NULL, tableHeader, "EXTTYPE"); if (!exttype) psAbort("cannot read table type"); Index: trunk/psModules/src/objects/pmSourceIO_MatchedRefs.c =================================================================== --- trunk/psModules/src/objects/pmSourceIO_MatchedRefs.c (revision 26190) +++ trunk/psModules/src/objects/pmSourceIO_MatchedRefs.c (revision 26260) @@ -120,4 +120,6 @@ psMetadataAdd (row, PS_LIST_TAIL, "X_CHIP", PS_DATA_F32, "x coord on chip", raw->chip->x); psMetadataAdd (row, PS_LIST_TAIL, "Y_CHIP", PS_DATA_F32, "y coord on chip", raw->chip->y); + psMetadataAdd (row, PS_LIST_TAIL, "X_CHIP_FIT",PS_DATA_F32, "x fitted coord on chip", ref->chip->x); + psMetadataAdd (row, PS_LIST_TAIL, "Y_CHIP_FIT",PS_DATA_F32, "y fitted coord on chip", ref->chip->y); psMetadataAdd (row, PS_LIST_TAIL, "X_FPA", PS_DATA_F32, "x coord on focal plane", raw->FP->x); psMetadataAdd (row, PS_LIST_TAIL, "Y_FPA", PS_DATA_F32, "y coord on focal plane", raw->FP->y); Index: trunk/psModules/src/objects/pmSourceVisual.c =================================================================== --- trunk/psModules/src/objects/pmSourceVisual.c (revision 26190) +++ trunk/psModules/src/objects/pmSourceVisual.c (revision 26260) @@ -73,4 +73,167 @@ float range; + range = graphdata.xmax - graphdata.xmin; + graphdata.xmax += 0.05*range; + graphdata.xmin -= 0.05*range; + range = graphdata.ymax - graphdata.ymin; + graphdata.ymax += 0.05*range; + graphdata.ymin -= 0.05*range; + + // better choice for range? + // graphdata.xmin = -17.0; + // graphdata.xmax = -9.0; + graphdata.ymin = -0.51; + graphdata.ymax = +0.51; + + KapaSetLimits (kapa1, &graphdata); + + KapaSetFont (kapa1, "helvetica", 14); + KapaBox (kapa1, &graphdata); + KapaSendLabel (kapa1, "PSF Mag", KAPA_LABEL_XM); + KapaSendLabel (kapa1, "Ap Mag - PSF Mag", KAPA_LABEL_YM); + + graphdata.color = KapaColorByName ("black"); + graphdata.ptype = 2; + graphdata.size = 0.5; + graphdata.style = 2; + graphdata.etype |= 0x01; + + KapaPrepPlot (kapa1, n, &graphdata); + KapaPlotVector (kapa1, n, x->data.F32, "x"); + KapaPlotVector (kapa1, n, y->data.F32, "y"); + KapaPlotVector (kapa1, n, dy->data.F32, "dym"); + KapaPlotVector (kapa1, n, dy->data.F32, "dyp"); + + psFree (x); + psFree (y); + psFree (dy); + + pmVisualAskUser(NULL); + return true; +} + +bool pmSourceVisualPlotPSFMetricSubpix (pmPSFtry *psfTry) { + + KapaSection section; // put the positive profile in one and the residuals in another? + Graphdata graphdata; + + if (!pmVisualIsVisual()) return true; + + if (kapa1 == -1) { + kapa1 = KapaOpenNamedSocket ("kapa", "pmSource:plots"); + if (kapa1 == -1) { + fprintf (stderr, "failure to open kapa; visual mode disabled\n"); + pmVisualSetVisual(false); + return false; + } + } + + KapaClearSections (kapa1); + KapaInitGraph (&graphdata); + + int n; + float range; + psVector *x = psVectorAllocEmpty (psfTry->sources->n, PS_TYPE_F32); + psVector *y = psVectorAllocEmpty (psfTry->sources->n, PS_TYPE_F32); + psVector *dy = psVectorAllocEmpty(psfTry->sources->n, PS_TYPE_F32); + + // section a: fractional-x pixel + section.dx = 1.0; + section.dy = 0.5; + section.x = 0.0; + section.y = 0.0; + section.name = NULL; + psStringAppend (§ion.name, "a1"); + KapaSetSection (kapa1, §ion); + psFree (section.name); + + graphdata.xmin = +32.0; + graphdata.xmax = -32.0; + graphdata.ymin = +32.0; + graphdata.ymax = -32.0; + + // construct the plot vectors + n = 0; + for (int i = 0; i < psfTry->sources->n; i++) { + if (psfTry->mask->data.PS_TYPE_VECTOR_MASK_DATA[i] & PSFTRY_MASK_ALL) continue; + + pmSource *source = psfTry->sources->data[i]; + x->data.F32[n] = source->modelEXT->params->data.F32[PM_PAR_XPOS] - (int)source->modelEXT->params->data.F32[PM_PAR_XPOS]; + + y->data.F32[n] = psfTry->metric->data.F32[i]; + dy->data.F32[n] = psfTry->metricErr->data.F32[i]; + graphdata.xmin = PS_MIN(graphdata.xmin, x->data.F32[n]); + graphdata.xmax = PS_MAX(graphdata.xmax, x->data.F32[n]); + graphdata.ymin = PS_MIN(graphdata.ymin, y->data.F32[n]); + graphdata.ymax = PS_MAX(graphdata.ymax, y->data.F32[n]); + n++; + } + x->n = y->n = dy->n = n; + + range = graphdata.xmax - graphdata.xmin; + graphdata.xmax += 0.05*range; + graphdata.xmin -= 0.05*range; + range = graphdata.ymax - graphdata.ymin; + graphdata.ymax += 0.05*range; + graphdata.ymin -= 0.05*range; + + // better choice for range? + // graphdata.xmin = -17.0; + // graphdata.xmax = -9.0; + graphdata.ymin = -0.51; + graphdata.ymax = +0.51; + + KapaSetLimits (kapa1, &graphdata); + + KapaSetFont (kapa1, "helvetica", 14); + KapaBox (kapa1, &graphdata); + KapaSendLabel (kapa1, "PSF Mag", KAPA_LABEL_XM); + KapaSendLabel (kapa1, "Ap Mag - PSF Mag", KAPA_LABEL_YM); + + graphdata.color = KapaColorByName ("black"); + graphdata.ptype = 2; + graphdata.size = 0.5; + graphdata.style = 2; + graphdata.etype |= 0x01; + + KapaPrepPlot (kapa1, n, &graphdata); + KapaPlotVector (kapa1, n, x->data.F32, "x"); + KapaPlotVector (kapa1, n, y->data.F32, "y"); + KapaPlotVector (kapa1, n, dy->data.F32, "dym"); + KapaPlotVector (kapa1, n, dy->data.F32, "dyp"); + + // *** section b: fractional-x pixel + section.dx = 1.0; + section.dy = 0.5; + section.x = 0.0; + section.y = 0.5; + section.name = NULL; + psStringAppend (§ion.name, "a2"); + KapaSetSection (kapa1, §ion); + psFree (section.name); + + graphdata.xmin = +32.0; + graphdata.xmax = -32.0; + graphdata.ymin = +32.0; + graphdata.ymax = -32.0; + + // construct the plot vectors + n = 0; + for (int i = 0; i < psfTry->sources->n; i++) { + if (psfTry->mask->data.PS_TYPE_VECTOR_MASK_DATA[i] & PSFTRY_MASK_ALL) continue; + + pmSource *source = psfTry->sources->data[i]; + x->data.F32[n] = source->modelEXT->params->data.F32[PM_PAR_YPOS] - (int)source->modelEXT->params->data.F32[PM_PAR_YPOS]; + + y->data.F32[n] = psfTry->metric->data.F32[i]; + dy->data.F32[n] = psfTry->metricErr->data.F32[i]; + graphdata.xmin = PS_MIN(graphdata.xmin, x->data.F32[n]); + graphdata.xmax = PS_MAX(graphdata.xmax, x->data.F32[n]); + graphdata.ymin = PS_MIN(graphdata.ymin, y->data.F32[n]); + graphdata.ymax = PS_MAX(graphdata.ymax, y->data.F32[n]); + n++; + } + x->n = y->n = dy->n = n; + range = graphdata.xmax - graphdata.xmin; graphdata.xmax += 0.05*range; Index: trunk/psModules/src/objects/pmSourceVisual.h =================================================================== --- trunk/psModules/src/objects/pmSourceVisual.h (revision 26190) +++ trunk/psModules/src/objects/pmSourceVisual.h (revision 26260) @@ -19,4 +19,6 @@ bool pmSourceVisualPSFModelResid (pmTrend2D *trend, psVector *x, psVector *y, psVector *param, psVector *mask); bool pmSourceVisualPlotPSFMetric (pmPSFtry *try); +bool pmSourceVisualPlotPSFMetricSubpix (pmPSFtry *try); + bool pmSourceVisualShowModelFit (pmSource *source); bool pmSourceVisualShowModelFits (pmPSF *psf, psArray *sources, psImageMaskType maskVal);