Index: /branches/eam_rel8_b2/psModules/src/objects/pmObjects.c =================================================================== --- /branches/eam_rel8_b2/psModules/src/objects/pmObjects.c (revision 5658) +++ /branches/eam_rel8_b2/psModules/src/objects/pmObjects.c (revision 5659) @@ -6,6 +6,6 @@ * @author EAM, IfA: significant modifications. * - * @version $Revision: 1.2.6.1 $ $Name: not supported by cvs2svn $ - * @date $Date: 2005-11-26 02:48:02 $ + * @version $Revision: 1.2.6.1.2.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2005-12-02 21:27:57 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -81,4 +81,7 @@ tmp->chisq = 0.0; tmp->nIter = 0; + tmp->radius = 0; + tmp->status = PM_MODEL_UNTRIED; + psS32 Nparams = pmModelParameterCount(type); if (Nparams == 0) { @@ -916,8 +919,16 @@ psImage *splane = NULL; int binX, binY; + bool status; + + psF32 SX_MAX = psMetadataLookupF32 (&status, metadata, "MOMENTS_SX_MAX"); + if (!status) + SX_MAX = 10.0; + psF32 SY_MAX = psMetadataLookupF32 (&status, metadata, "MOMENTS_SY_MAX"); + if (!status) + SY_MAX = 10.0; // construct a sigma-plane image // psImageAlloc does zero the data - splane = psImageAlloc (NPIX/SCALE, NPIX/SCALE, PS_TYPE_F32); + splane = psImageAlloc (SX_MAX/SCALE, SY_MAX/SCALE, PS_TYPE_F32); for (int i = 0; i < splane->numRows; i++) { @@ -1031,6 +1042,6 @@ psfClump.dY = stats->clippedStdev; - psTrace (".pmObjects.pmSourceRoughClass", 2, "clump X, Y: %f, %f\n", psfClump.X, psfClump.Y); - psTrace (".pmObjects.pmSourceRoughClass", 2, "clump DX, DY: %f, %f\n", psfClump.dX, psfClump.dY); + psLogMsg ("pmObjects", 4, "clump X, Y: %f, %f\n", psfClump.X, psfClump.Y); + psLogMsg ("pmObjects", 4, "clump DX, DY: %f, %f\n", psfClump.dX, psfClump.dY); // these values should be pushed on the metadata somewhere @@ -1160,13 +1171,54 @@ } - psTrace (".pmObjects.pmSourceRoughClass", 2, "Nstar: %3d\n", Nstar); - psTrace (".pmObjects.pmSourceRoughClass", 2, "Npsf: %3d\n", Npsf); - psTrace (".pmObjects.pmSourceRoughClass", 2, "Ngal: %3d\n", Ngal); - psTrace (".pmObjects.pmSourceRoughClass", 2, "Nsatstar: %3d\n", Nsatstar); - psTrace (".pmObjects.pmSourceRoughClass", 2, "Nsat: %3d\n", Nsat); - psTrace (".pmObjects.pmSourceRoughClass", 2, "Ncr: %3d\n", Ncr); + psLogMsg ("pmObjects", 4, "Nstar: %3d\n", Nstar); + psLogMsg ("pmObjects", 4, "Npsf: %3d\n", Npsf); + psLogMsg ("pmObjects", 4, "Ngal: %3d\n", Ngal); + psLogMsg ("pmObjects", 4, "Nsatstar: %3d\n", Nsatstar); + psLogMsg ("pmObjects", 4, "Nsat: %3d\n", Nsat); + psLogMsg ("pmObjects", 4, "Ncr: %3d\n", Ncr); + + // optional dump of all rough source data + { + bool status; + char *output = psMetadataLookupPtr (&status, metadata, "MOMENTS_OUTPUT_FILE"); + if (status && (output != NULL) && (output[0])) + { + pmMomentsWriteText (sources, output); + psFree (output); + } + } psTrace(__func__, 3, "---- %s(%d) end ----\n", __func__, rc); return(rc); +} + +// write the moments to an output file +bool pmMomentsWriteText (psArray *sources, char *filename) +{ + + int i; + FILE *f; + pmSource *source; + + f = fopen (filename, "w"); + if (f == NULL) { + psLogMsg ("pmMomentsWriteText", 3, "can't open output file for moments%s\n", filename); + return false; + } + + for (i = 0; i < sources->n; i++) { + source = (pmSource *) sources->data[i]; + if (source->moments == NULL) + continue; + fprintf (f, "%5d %5d %7.1f %7.1f %7.1f %6.3f %6.3f %8.1f %7.1f %7.1f %7.1f %4d %2d\n", + source->peak->x, source->peak->y, source->peak->counts, + source->moments->x, source->moments->y, + source->moments->Sx, source->moments->Sy, + source->moments->Sum, source->moments->Peak, + source->moments->Sky, source->moments->SN, + source->moments->nPixels, source->type); + } + fclose (f); + return true; } @@ -1575,5 +1627,4 @@ psBool fitStatus = true; psBool onPic = true; - psBool rc = true; // XXX EAM : is it necessary for the mask & weight to exist? the @@ -1602,4 +1653,5 @@ psTrace (".pmObjects.pmSourceFitModel", 4, "insufficient valid pixels\n"); psTrace(__func__, 3, "---- %s(false) end ----\n", __func__); + model->status = PM_MODEL_BADARGS; return(false); } @@ -1642,5 +1694,5 @@ } - // XXX EAM : covar must be F64? + // covar must be F64 psImage *covar = psImageAlloc (params->n, params->n, PS_TYPE_F64); @@ -1653,19 +1705,10 @@ } - // XXX EAM: we need to do something (give an error?) if rc is false - // XXX EAM: psMinimizeLMChi2 does not check convergence - - // XXX models can go insane: reject these - onPic &= (params->data.F32[2] >= source->pixels->col0); - onPic &= (params->data.F32[2] < source->pixels->col0 + source->pixels->numCols); - onPic &= (params->data.F32[3] >= source->pixels->row0); - onPic &= (params->data.F32[3] < source->pixels->row0 + source->pixels->numRows); - - // XXX EAM: save the resulting chisq, nDOF, nIter + // save the resulting chisq, nDOF, nIter model->chisq = myMin->value; model->nIter = myMin->iter; model->nDOF = y->n - nParams; - // XXX EAM get the Gauss-Newton distance for fixed model parameters + // get the Gauss-Newton distance for fixed model parameters if (paramMask != NULL) { psVector *delta = psVectorAlloc (params->n, PS_TYPE_F64); @@ -1679,4 +1722,20 @@ } + // set the model success or failure status + if (!fitStatus) { + model->status = PM_MODEL_NONCONVERGE; + } else { + model->status = PM_MODEL_SUCCESS; + } + + // models can go insane: save this status instead + onPic &= (params->data.F32[2] >= source->pixels->col0); + onPic &= (params->data.F32[2] < source->pixels->col0 + source->pixels->numCols); + onPic &= (params->data.F32[3] >= source->pixels->row0); + onPic &= (params->data.F32[3] < source->pixels->row0 + source->pixels->numRows); + if (!onPic) { + model->status = PM_MODEL_OFFIMAGE; + } + psFree(x); psFree(y); @@ -1686,7 +1745,6 @@ psFree(paramMask); - rc = (onPic && fitStatus); psTrace(__func__, 3, "---- %s(%d) end ----\n", __func__, rc); - return(rc); + return(onPic && fitStatus); }