Index: branches/eam_branches/ipp-20110213/psModules/src/objects/pmModel.c =================================================================== --- branches/eam_branches/ipp-20110213/psModules/src/objects/pmModel.c (revision 30775) +++ branches/eam_branches/ipp-20110213/psModules/src/objects/pmModel.c (revision 30784) @@ -67,4 +67,5 @@ tmp->chisqNorm = NAN; tmp->nDOF = 0; + tmp->nPar = 0; tmp->nPix = 0; tmp->nIter = 0; Index: branches/eam_branches/ipp-20110213/psModules/src/objects/pmModel.h =================================================================== --- branches/eam_branches/ipp-20110213/psModules/src/objects/pmModel.h (revision 30775) +++ branches/eam_branches/ipp-20110213/psModules/src/objects/pmModel.h (revision 30784) @@ -39,5 +39,6 @@ float magErr; ///< integrated model magnitude error int nPix; ///< number of pixels used for fit - int nDOF; ///< number of degrees of freedom + int nPar; ///< number of parameters in fit + int nDOF; ///< number of degrees of freedom (nDOF = nPix - nPar) int nIter; ///< number of iterations to reach min pmModelStatus flags; ///< model status flags Index: branches/eam_branches/ipp-20110213/psModules/src/objects/pmPCMdata.c =================================================================== --- branches/eam_branches/ipp-20110213/psModules/src/objects/pmPCMdata.c (revision 30775) +++ branches/eam_branches/ipp-20110213/psModules/src/objects/pmPCMdata.c (revision 30784) @@ -254,5 +254,6 @@ pcm->nPix = nPix; - pcm->nDOF = nPix - nParams - 1; + pcm->nPar = nParams; + pcm->nDOF = nPix - nParams; return pcm; @@ -341,4 +342,5 @@ return false; } + pcm->nPar = nParams; pcm->nDOF = pcm->nPix - nParams; Index: branches/eam_branches/ipp-20110213/psModules/src/objects/pmPCMdata.h =================================================================== --- branches/eam_branches/ipp-20110213/psModules/src/objects/pmPCMdata.h (revision 30775) +++ branches/eam_branches/ipp-20110213/psModules/src/objects/pmPCMdata.h (revision 30784) @@ -33,4 +33,5 @@ psMinConstraint *constraint; int nPix; + int nPar; int nDOF; } pmPCMdata; Index: branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceFitModel.c =================================================================== --- branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceFitModel.c (revision 30775) +++ branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceFitModel.c (revision 30784) @@ -236,16 +236,19 @@ model->covar = psMemIncrRefCounter(covar); } + model->nIter = myMin->iter; + model->nPar = nParams; + psTrace ("psModules.objects", 4, "niter: %d, chisq: %f", myMin->iter, myMin->value); // save the resulting chisq, nDOF, nIter + // NOTE: if (!options->poissonErrors) chisq will be wrong : recalculate if (options->poissonErrors) { model->chisq = myMin->value; model->nPix = y->n; - model->nDOF = y->n - nParams; + model->nDOF = y->n - model->nPar; model->chisqNorm = model->chisq / model->nDOF; } else { - pmSourceChisq (model, source->pixels, source->maskObj, source->variance, maskVal, options->covarFactor, nParams); - } - model->nIter = myMin->iter; + pmSourceChisqUnsubtracted (source, model, maskVal); + } // set the model success or failure status Index: branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceFitPCM.c =================================================================== --- branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceFitPCM.c (revision 30775) +++ branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceFitPCM.c (revision 30784) @@ -84,4 +84,6 @@ } } + pcm->modelConv->nIter = myMin->iter; + pcm->modelConv->nPar = pcm->nPar; // save the resulting chisq, nDOF, nIter @@ -92,7 +94,7 @@ pcm->modelConv->chisqNorm = pcm->modelConv->chisq / pcm->modelConv->nDOF; } else { - pmSourceChisq (pcm->modelConv, source->pixels, source->maskObj, source->variance, maskVal, fitOptions->covarFactor, pcm->nPix - pcm->nDOF - 1); + // xxx this is wrong because it does not convolve with the psf + pmSourceChisqUnsubtracted (source, pcm->modelConv, maskVal); } - pcm->modelConv->nIter = myMin->iter; // set the model success or failure status Index: branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceFitSet.c =================================================================== --- branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceFitSet.c (revision 30775) +++ branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceFitSet.c (revision 30784) @@ -335,4 +335,7 @@ psTrace ("psModules.objects", 4, " src %d", i); + model->nIter = myMin->iter; + // model->nPar is set by pmSourceFitSetMasks + // save the resulting chisq, nDOF, nIter // these are not unique for any one source @@ -340,10 +343,9 @@ model->chisq = myMin->value; model->nPix = nPix; - model->nDOF = nPix - model->params->n; + model->nDOF = nPix - model->nPar; model->chisqNorm = model->chisq / model->nDOF; } else { - pmSourceChisq (model, source->pixels, source->maskObj, source->variance, maskVal, options->covarFactor, model->params->n); + pmSourceChisqUnsubtracted (source, model, maskVal); } - model->nIter = myMin->iter; // set the model success or failure status @@ -399,4 +401,5 @@ for (int i = 0; i < set->paramSet->n; i++) { psVector *paramOne = set->paramSet->data[i]; + pmModel *modelOne = set->modelSet->data[i]; switch (mode) { @@ -406,4 +409,5 @@ if (j == PM_PAR_I0) continue; constraint->paramMask->data.PS_TYPE_VECTOR_MASK_DATA[n + j] = 1; + modelOne->nPar = 1; } break; @@ -415,4 +419,5 @@ if (j == PM_PAR_I0) continue; constraint->paramMask->data.PS_TYPE_VECTOR_MASK_DATA[n + j] = 1; + modelOne->nPar = 3; } break; @@ -420,4 +425,5 @@ // EXT model fits all params (except sky) constraint->paramMask->data.PS_TYPE_VECTOR_MASK_DATA[n + PM_PAR_SKY] = 1; + modelOne->nPar = paramOne->n - 1; break; default: Index: branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceOutputs.c =================================================================== --- branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceOutputs.c (revision 30775) +++ branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourceOutputs.c (revision 30784) @@ -72,5 +72,5 @@ } - *nImageOverlap = 1; + *nImageOverlap = psMetadataLookupS32 (&status2, header, "NINPUTS"); return true; Index: branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.c =================================================================== --- branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.c (revision 30775) +++ branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.c (revision 30784) @@ -741,6 +741,6 @@ # endif -// determine chisq, etc for linear normalization-only fit -bool pmSourceChisq (pmModel *model, psImage *image, psImage *mask, psImage *variance, psImageMaskType maskVal, const float covarFactor, int nParams) +// determine chisq, nPix, nDOF, chisqNorm : model->nPar must be set +bool pmSourceChisq (pmModel *model, psImage *image, psImage *mask, psImage *variance, psImageMaskType maskVal) { PS_ASSERT_PTR_NON_NULL(model, false); @@ -757,12 +757,10 @@ if (variance->data.F32[j][i] <= 0) continue; - // dC += PS_SQR (image->data.F32[j][i]) / (covarFactor * variance->data.F32[j][i]); dC += PS_SQR (image->data.F32[j][i]) / variance->data.F32[j][i]; Npix ++; } } - model->nPix = Npix; - model->nDOF = Npix - nParams - 1; + model->nDOF = Npix - model->nPar; model->chisq = dC; model->chisqNorm = dC / model->nDOF; @@ -771,4 +769,54 @@ } + +// return source aperture magnitude +bool pmSourceChisqUnsubtracted (pmSource *source, pmModel *model, psImageMaskType maskVal) +{ + PS_ASSERT_PTR_NON_NULL(source, false); + PS_ASSERT_PTR_NON_NULL(model, false); + + float dC = 0.0; + int Npix = 0; + + // the model function returns the source flux at a position + psVector *coord = psVectorAlloc(2, PS_TYPE_F32); + + psVector *params = model->params; + psImage *image = source->pixels; + psImage *mask = source->maskObj; + psImage *variance = source->variance; + + int dX = image->col0; + int dY = image->row0; + + for (int iy = 0; iy < image->numRows; iy++) { + for (int ix = 0; ix < image->numCols; ix++) { + + // skip pixels which are masked + if (mask->data.PS_TYPE_IMAGE_MASK_DATA[iy][ix] & maskVal) continue; + + if (variance->data.F32[iy][ix] <= 0) continue; + + coord->data.F32[0] = (psF32) ix + dX + 0.5; + coord->data.F32[1] = (psF32) iy + dY + 0.5; + + // for the full model, add all points + float value = model->modelFunc (NULL, params, coord); + + // fprintf (stderr, "%d, %d : %f, %f : %f - %f : %f\n", + // ix, iy, coord->data.F32[0], coord->data.F32[1], image->data.F32[iy][ix], value, dC); + + dC += PS_SQR (image->data.F32[iy][ix] - value) / variance->data.F32[iy][ix]; + Npix ++; + } + } + model->nPix = Npix; + model->nDOF = Npix - model->nPar; + model->chisq = dC; + model->chisqNorm = dC / model->nDOF; + + psFree (coord); + return (true); +} double pmSourceModelWeight(const pmSource *Mi, int term, const bool unweighted_sum, const float covarFactor, psImageMaskType maskVal) Index: branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.h =================================================================== --- branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.h (revision 30775) +++ branches/eam_branches/ipp-20110213/psModules/src/objects/pmSourcePhotometry.h (revision 30784) @@ -69,5 +69,6 @@ bool pmSourcePixelWeight (pmSource *source, pmModel *model, psImage *mask, psImageMaskType maskVal, float radius); -bool pmSourceChisq (pmModel *model, psImage *image, psImage *mask, psImage *weight, psImageMaskType maskVal, const float covarFactor, int nParams); +bool pmSourceChisq (pmModel *model, psImage *image, psImage *mask, psImage *weight, psImageMaskType maskVal); +bool pmSourceChisqUnsubtracted (pmSource *source, pmModel *model, psImageMaskType maskVal); bool pmSourceMeasureDiffStats (pmSource *source, psImageMaskType maskVal, psImageMaskType markVal);