Index: /branches/eam_branch_20070921/psModules/src/camera/pmFPAfileIO.c =================================================================== --- /branches/eam_branch_20070921/psModules/src/camera/pmFPAfileIO.c (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/camera/pmFPAfileIO.c (revision 14989) @@ -774,4 +774,5 @@ } pmFPAview *thisView = pmFPAAddSourceFromHeader (file->fpa, phu, file->format); + assert (thisView); // XXX we are having some trouble with input psf files not having the Cell and fpa names matching. psFree (thisView); psFree (phu); Index: /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_GAUSS.c =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_GAUSS.c (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_GAUSS.c (revision 14989) @@ -274,6 +274,6 @@ out[i] = in[i]; } else { - psPolynomial2D *poly = psf->params->data[i]; - out[i] = psPolynomial2DEval(poly, in[PM_PAR_XPOS], in[PM_PAR_YPOS]); + pmTrend2D *trend = psf->params->data[i]; + out[i] = pmTrend2DEval(trend, in[PM_PAR_XPOS], in[PM_PAR_YPOS]); } } @@ -331,7 +331,6 @@ if (i == PM_PAR_XPOS) continue; if (i == PM_PAR_YPOS) continue; - psPolynomial2D *poly = psf->params->data[i]; - assert (poly); - PAR[i] = psPolynomial2DEval(poly, Xo, Yo); + pmTrend2D *trend = psf->params->data[i]; + PAR[i] = pmTrend2DEval(trend, Xo, Yo); } Index: /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_PGAUSS.c =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_PGAUSS.c (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_PGAUSS.c (revision 14989) @@ -328,6 +328,6 @@ out[i] = in[i]; } else { - psPolynomial2D *poly = psf->params->data[i]; - out[i] = psPolynomial2DEval(poly, in[PM_PAR_XPOS], in[PM_PAR_YPOS]); + pmTrend2D *trend = psf->params->data[i]; + out[i] = pmTrend2DEval(trend, in[PM_PAR_XPOS], in[PM_PAR_YPOS]); } } @@ -384,7 +384,6 @@ if (i == PM_PAR_XPOS) continue; if (i == PM_PAR_YPOS) continue; - psPolynomial2D *poly = psf->params->data[i]; - assert (poly); - PAR[i] = psPolynomial2DEval(poly, Xo, Yo); + pmTrend2D *trend = psf->params->data[i]; + PAR[i] = pmTrend2DEval(trend, Xo, Yo); } Index: /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_QGAUSS.c =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_QGAUSS.c (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_QGAUSS.c (revision 14989) @@ -358,6 +358,6 @@ out[i] = in[i]; } else { - psPolynomial2D *poly = psf->params->data[i]; - out[i] = psPolynomial2DEval(poly, in[PM_PAR_XPOS], in[PM_PAR_YPOS]); + pmTrend2D *trend = psf->params->data[i]; + out[i] = pmTrend2DEval(trend, in[PM_PAR_XPOS], in[PM_PAR_YPOS]); } } @@ -412,7 +412,6 @@ if (i == PM_PAR_XPOS) continue; if (i == PM_PAR_YPOS) continue; - psPolynomial2D *poly = psf->params->data[i]; - assert (poly); - PAR[i] = psPolynomial2DEval(poly, Xo, Yo); + pmTrend2D *trend = psf->params->data[i]; + PAR[i] = pmTrend2DEval(trend, Xo, Yo); } Index: /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_RGAUSS.c =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_RGAUSS.c (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_RGAUSS.c (revision 14989) @@ -351,6 +351,6 @@ out[i] = in[i]; } else { - psPolynomial2D *poly = psf->params->data[i]; - out[i] = psPolynomial2DEval(poly, in[PM_PAR_XPOS], in[PM_PAR_YPOS]); + pmTrend2D *trend = psf->params->data[i]; + out[i] = pmTrend2DEval(trend, in[PM_PAR_XPOS], in[PM_PAR_YPOS]); } } @@ -405,7 +405,6 @@ if (i == PM_PAR_XPOS) continue; if (i == PM_PAR_YPOS) continue; - psPolynomial2D *poly = psf->params->data[i]; - assert (poly); - PAR[i] = psPolynomial2DEval(poly, Xo, Yo); + pmTrend2D *trend = psf->params->data[i]; + PAR[i] = pmTrend2DEval(trend, Xo, Yo); } Index: /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_SERSIC.c =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_SERSIC.c (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/models/pmModel_SERSIC.c (revision 14989) @@ -343,6 +343,6 @@ out[i] = in[i]; } else { - psPolynomial2D *poly = psf->params->data[i]; - out[i] = psPolynomial2DEval(poly, in[PM_PAR_XPOS], in[PM_PAR_YPOS]); + pmTrend2D *trend = psf->params->data[i]; + out[i] = pmTrend2DEval(trend, in[PM_PAR_XPOS], in[PM_PAR_YPOS]); } } @@ -397,7 +397,6 @@ if (i == PM_PAR_XPOS) continue; if (i == PM_PAR_YPOS) continue; - psPolynomial2D *poly = psf->params->data[i]; - assert (poly); - PAR[i] = psPolynomial2DEval(poly, Xo, Yo); + pmTrend2D *trend = psf->params->data[i]; + PAR[i] = pmTrend2DEval(trend, Xo, Yo); } Index: /branches/eam_branch_20070921/psModules/src/objects/pmPSF.c =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/pmPSF.c (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/pmPSF.c (revision 14989) @@ -6,6 +6,6 @@ * @author EAM, IfA * - * @version $Revision: 1.28.2.4 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-09-22 13:47:28 $ + * @version $Revision: 1.28.2.5 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-09-24 03:00:44 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -48,5 +48,4 @@ if (!options) return; - psFree (options->image); psFree (options->stats); return; @@ -60,5 +59,4 @@ options->type = 0; - options->image = NULL; options->stats = NULL; @@ -66,4 +64,8 @@ options->psfTrendNx = 0; options->psfTrendNy = 0; + options->psfFieldNx = 0; + options->psfFieldNy = 0; + options->psfFieldXo = 0; + options->psfFieldYo = 0; options->poissonErrorsPhotLMM = true; @@ -84,4 +86,5 @@ psFree (psf->ChiTrend); + psFree (psf->psfTrendStats); psFree (psf->ApTrend); psFree (psf->FluxScale); @@ -164,6 +167,15 @@ // PM_PAR_XPOS, etc) - psf->Nx = options->image->numCols; - psf->Ny = options->image->numRows; + psImageBinning *binning = psImageBinningAlloc(); + binning->nXruff = options->psfTrendNx; + binning->nYruff = options->psfTrendNy; + binning->nXfine = options->psfFieldNx; + binning->nYfine = options->psfFieldNy; + psImageBinningSetSkipByOffset (binning, options->psfFieldXo, options->psfFieldYo); + + psf->fieldNx = options->psfFieldNx; + psf->fieldNy = options->psfFieldNy; + psf->fieldXo = options->psfFieldXo; + psf->fieldYo = options->psfFieldYo; // define the parameter trends @@ -175,7 +187,8 @@ if (i == PM_PAR_YPOS) continue; - psf->params->data[i] = pmTrend2DAlloc (options->psfTrendMode, options->image, options->psfTrendNx, options->psfTrendNx, options->stats); + psf->params->data[i] = pmTrend2DNoImageAlloc (options->psfTrendMode, binning, options->stats); } } + psFree (binning); return psf; } Index: /branches/eam_branch_20070921/psModules/src/objects/pmPSF.h =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/pmPSF.h (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/pmPSF.h (revision 14989) @@ -6,6 +6,6 @@ * @author EAM, IfA * - * @version $Revision: 1.16.2.3 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-09-22 13:47:28 $ + * @version $Revision: 1.16.2.4 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-09-24 03:00:44 $ * Copyright 2004 Maui High Performance Computing Center, University of Hawaii */ @@ -47,5 +47,8 @@ int nPSFstars; ///< number of stars used to measure PSF int nApResid; ///< number of stars used to measure ApResid - int Nx, Ny; ///< reference image dimensions + int fieldNx; + int fieldNy; + int fieldXo; + int fieldYo; bool poissonErrorsPhotLMM; ///< use poission errors for non-linear model fitting bool poissonErrorsPhotLin; ///< use poission errors for linear model fitting @@ -58,14 +61,16 @@ typedef struct { pmModelType type; - psImage *image; // image for which the PSF is defined psStats *stats; // psfTrend clipping stats pmTrend2DMode psfTrendMode; int psfTrendNx; int psfTrendNy; + int psfFieldNx; + int psfFieldNy; + int psfFieldXo; + int psfFieldYo; bool poissonErrorsPhotLMM; ///< use poission errors for non-linear model fitting bool poissonErrorsPhotLin; ///< use poission errors for linear model fitting bool poissonErrorsParams; ///< use poission errors for model parameter fitting float radius; - } pmPSFOptions; Index: /branches/eam_branch_20070921/psModules/src/objects/pmPSF_IO.c =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/pmPSF_IO.c (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/pmPSF_IO.c (revision 14989) @@ -6,6 +6,6 @@ * @author EAM, IfA * - * @version $Revision: 1.22.2.3 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-09-22 13:47:28 $ + * @version $Revision: 1.22.2.4 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-09-24 03:00:44 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -163,4 +163,11 @@ if (!analysis) return false; + // select the current recipe + psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, "PSPHOT"); + if (!recipe) { + psError(PS_ERR_UNKNOWN, false, "missing recipe %s", "PSPHOT"); + return false; + } + // write a PHU? (only if input image is MEF) // write a header? (only if this is the first readout for cell) @@ -260,6 +267,21 @@ psMetadataAdd (header, PS_LIST_TAIL, "PSF_NPAR", PS_DATA_S32, "PSF model parameter count", nPar); - psMetadataAddS32 (header, PS_LIST_TAIL, "IMAXIS1", 0, "Image X Size", psf->Nx); - psMetadataAddS32 (header, PS_LIST_TAIL, "IMAXIS2", 0, "Image Y Size", psf->Ny); + psMetadataAddS32 (header, PS_LIST_TAIL, "IMAXIS1", 0, "Image X Size", psf->fieldNx); + psMetadataAddS32 (header, PS_LIST_TAIL, "IMAXIS2", 0, "Image Y Size", psf->fieldNy); + psMetadataAddS32 (header, PS_LIST_TAIL, "IMREF1", 0, "Image X Ref", psf->fieldXo); + psMetadataAddS32 (header, PS_LIST_TAIL, "IMREF2", 0, "Image Y Ref", psf->fieldYo); + + // extract PSF Clump info + pmPSFClump psfClump; + + psfClump.X = psMetadataLookupF32 (&status, recipe, "PSF.CLUMP.X"); assert (status); + psfClump.Y = psMetadataLookupF32 (&status, recipe, "PSF.CLUMP.Y"); assert (status); + psfClump.dX = psMetadataLookupF32 (&status, recipe, "PSF.CLUMP.DX"); assert (status); + psfClump.dY = psMetadataLookupF32 (&status, recipe, "PSF.CLUMP.DY"); assert (status); + + psMetadataAddF32 (header, PS_LIST_TAIL, "PSF_CLX", 0, "psf clump center", psfClump.X); + psMetadataAddF32 (header, PS_LIST_TAIL, "PSF_CLY", 0, "psf clump center", psfClump.Y); + psMetadataAddF32 (header, PS_LIST_TAIL, "PSF_CLDX", 0, "psf clump size", psfClump.dX); + psMetadataAddF32 (header, PS_LIST_TAIL, "PSF_CLDY", 0, "psf clump size", psfClump.dY); // save the dimensions of each parameter @@ -306,6 +328,6 @@ // write the image components into a table: this is needed because they may each be a different size psImageMap *map = trend->map; - for (int ix = 0; ix <= map->map->numCols; ix++) { - for (int iy = 0; iy <= map->map->numRows; iy++) { + for (int ix = 0; ix < map->map->numCols; ix++) { + for (int iy = 0; iy < map->map->numRows; iy++) { psMetadata *row = psMetadataAlloc (); psMetadataAddS32 (row, PS_LIST_TAIL, "MODEL_TERM", 0, "", i); @@ -313,5 +335,5 @@ psMetadataAddS32 (row, PS_LIST_TAIL, "Y_POWER", 0, "", iy); psMetadataAddF32 (row, PS_LIST_TAIL, "VALUE", 0, "", map->map->data.F32[iy][ix]); - psMetadataAddF32 (row, PS_LIST_TAIL, "ERROR", 0, "", map->error->data.F32[ix][iy]); + psMetadataAddF32 (row, PS_LIST_TAIL, "ERROR", 0, "", map->error->data.F32[iy][ix]); psMetadataAddU8 (row, PS_LIST_TAIL, "MASK", 0, "", 0); // no cells are masked @@ -518,4 +540,11 @@ psTrace ("psModules.objects", 5, "read psf model for %s\n", file->filename); + // select the current recipe + psMetadata *recipe = psMetadataLookupPtr (NULL, config->recipes, "PSPHOT"); + if (!recipe) { + psError(PS_ERR_UNKNOWN, false, "missing recipe %s", "PSPHOT"); + return false; + } + // Menu of EXTNAME rules psMetadata *menu = psMetadataLookupMetadata(&status, file->camera, "EXTNAME.RULES"); @@ -560,11 +589,29 @@ } + // psf clump data + pmPSFClump psfClump; + + psfClump.X = psMetadataLookupF32 (&status, header, "PSF_CLX" ); assert(status); + psfClump.Y = psMetadataLookupF32 (&status, header, "PSF_CLY" ); assert(status); + psfClump.dX = psMetadataLookupF32 (&status, header, "PSF_CLDX"); assert(status); + psfClump.dY = psMetadataLookupF32 (&status, header, "PSF_CLDY"); assert(status); + + psMetadataAddF32 (recipe, PS_LIST_TAIL, "PSF.CLUMP.X" , 0, "psf clump center", psfClump.X); + psMetadataAddF32 (recipe, PS_LIST_TAIL, "PSF.CLUMP.Y" , 0, "psf clump center", psfClump.Y); + psMetadataAddF32 (recipe, PS_LIST_TAIL, "PSF.CLUMP.DX", 0, "psf clump size", psfClump.dX); + psMetadataAddF32 (recipe, PS_LIST_TAIL, "PSF.CLUMP.DY", 0, "psf clump size", psfClump.dY); + options->poissonErrorsPhotLMM = psMetadataLookupBool (&status, header, "ERR_LMM"); options->poissonErrorsPhotLin = psMetadataLookupBool (&status, header, "ERR_LIN"); options->poissonErrorsParams = psMetadataLookupBool (&status, header, "ERR_PAR"); - // XXX this is a bit ad-hoc: determine the PSF trend mode by looking at one of the parameters - // this could allow the psf model to use different trend modes for different parameters, but we - // would need to change the alloc somewhat + options->psfFieldNx = psMetadataLookupS32 (&status, header, "IMAXIS1"); + options->psfFieldNy = psMetadataLookupS32 (&status, header, "IMAXIS2"); + options->psfFieldXo = psMetadataLookupS32 (&status, header, "IMREF1"); + options->psfFieldYo = psMetadataLookupS32 (&status, header, "IMREF2"); + + psImageBinning *binning = psImageBinningAlloc(); + binning->nXfine = options->psfFieldNx; + binning->nYfine = options->psfFieldNy; // we determine the PSF parameter polynomials from the MD-defined polynomials @@ -576,20 +623,18 @@ psAbort("mismatch model par count"); - int nXimage = psMetadataLookupS32 (&status, header, "IMAXIS1"); - int nYimage = psMetadataLookupS32 (&status, header, "IMAXIS2"); - // XXX don't require this to exist? 0,0 if not defined? - // load the trend mode and dimensions of each parameter for (int i = 0; i < nPar; i++) { char name[9]; snprintf (name, 9, "PAR%02d_NX", i); - int nXorder = psMetadataLookupS32 (&status, header, name); + binning->nXruff = psMetadataLookupS32 (&status, header, name); if (!status) continue; // not all parameters are defined + snprintf (name, 9, "PAR%02d_NY", i); - int nYorder = psMetadataLookupS32 (&status, header, name); + binning->nYruff = psMetadataLookupS32 (&status, header, name); if (!status) { psError(PS_ERR_UNKNOWN, true, "inconsistent PSF header: NX defined for PAR %d, but not NY", i); return false; } + snprintf (name, 9, "PAR%02d_MD", i); char *modeName = psMetadataLookupStr (&status, header, name); @@ -599,8 +644,12 @@ } pmTrend2DMode psfTrendMode = pmTrend2DModeFromString (modeName); - // XXX default to POLY_ORD if not defined? - - psf->params->data[i] = pmTrend2DFieldAlloc (psfTrendMode, nXimage, nYimage, nXorder, nYorder, NULL); - } + if (psfTrendMode == PM_TREND_NONE) { + psfTrendMode = PM_TREND_POLY_ORD; + } + + psImageBinningSetSkipByOffset (binning, options->psfFieldXo, options->psfFieldYo); + psf->params->data[i] = pmTrend2DNoImageAlloc (psfTrendMode, binning, NULL); + } + psFree (binning); // other required information describing the PSF @@ -622,5 +671,4 @@ int xPow = psMetadataLookupS32 (&status, row, "X_POWER"); int yPow = psMetadataLookupS32 (&status, row, "Y_POWER"); - // XXX sanity check here pmTrend2D *trend = psf->params->data[iPar]; @@ -632,8 +680,20 @@ if (trend->mode == PM_TREND_MAP) { psImageMap *map = trend->map; + assert (map); + assert (map->map); + assert (map->error); + assert (xPow >= 0); + assert (yPow >= 0); + assert (xPow < map->map->numCols); + assert (yPow < map->map->numRows); map->map->data.F32[yPow][xPow] = psMetadataLookupF32 (&status, row, "VALUE"); map->error->data.F32[yPow][xPow] = psMetadataLookupF32 (&status, row, "ERROR"); } else { psPolynomial2D *poly = trend->poly; + assert (poly); + assert (xPow >= 0); + assert (yPow >= 0); + assert (xPow <= poly->nX); + assert (yPow <= poly->nY); poly->coeff[xPow][yPow] = psMetadataLookupF32 (&status, row, "VALUE"); poly->coeffErr[xPow][yPow] = psMetadataLookupF32 (&status, row, "ERROR"); Index: /branches/eam_branch_20070921/psModules/src/objects/pmSourcePhotometry.c =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/pmSourcePhotometry.c (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/pmSourcePhotometry.c (revision 14989) @@ -3,6 +3,6 @@ * @author EAM, IfA; GLG, MHPCC * - * @version $Revision: 1.31 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-09-21 02:46:46 $ + * @version $Revision: 1.31.2.1 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-09-24 03:00:44 $ * * Copyright 2004 Maui High Performance Computing Center, University of Hawaii @@ -81,6 +81,8 @@ // we must have a valid model model = pmSourceGetModel (&isPSF, source); - if (model == NULL) - return false; + if (model == NULL) { + psTrace ("psModules.objects", 3, "fail mag : no valid model"); + return false; + } if (model->dparams->data.F32[PM_PAR_I0] > 0) { @@ -109,12 +111,14 @@ // measure EXT model photometry status = pmSourcePhotometryModel (&source->extMag, source->modelEXT); + // for PSFs, correct both apMag and psfMag to same system, consistent with infinite flux star in aperture RADIUS - if ((mode & PM_SOURCE_PHOT_APCORR) && isPSF && psf) { - // convert to the equivalent 2D model? + if ((mode & PM_SOURCE_PHOT_APCORR) && isPSF && psf && psf->ApTrend) { source->psfMag += pmTrend2DEval (psf->ApTrend, x, y); } - if (!isfinite(SN) || (SN < AP_MIN_SN)) - return false; + if (!isfinite(SN) || (SN < AP_MIN_SN)) { + psTrace ("psModules.objects", 3, "fail mag : bad SN: %f (limit: %f)", SN, AP_MIN_SN); + return false; + } // replace source flux @@ -209,4 +213,5 @@ status = pmSourcePhotometryAper (&source->apMag, model, flux, source->maskObj, maskVal); if (!status) { + psTrace ("psModules.objects", 3, "fail mag : bad Ap Mag"); psErrorCode last = psErrorCodeLast(); if (last == PM_ERR_PHOTOM) { Index: /branches/eam_branch_20070921/psModules/src/objects/pmTrend2D.c =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/pmTrend2D.c (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/pmTrend2D.c (revision 14989) @@ -3,6 +3,6 @@ * @author EAM, IfA * - * @version $Revision: 1.2.2.3 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-09-22 13:47:28 $ + * @version $Revision: 1.2.2.4 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-09-24 03:00:44 $ * * Copyright 2004 Institute for Astronomy, University of Hawaii @@ -32,44 +32,5 @@ { assert (image); - assert (stats); - - pmTrend2D *trend = (pmTrend2D *) psAlloc(sizeof(pmTrend2D)); - psMemSetDeallocator(trend, (psFreeFunc) pmTrend2DFree); - - trend->map = NULL; - trend->poly = NULL; - trend->stats = psMemIncrRefCounter (stats); - trend->mode = mode; - - switch (mode) { - case PM_TREND_POLY_ORD: - trend->poly = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, nXtrend, nYtrend); - break; - - case PM_TREND_POLY_CHEB: - trend->poly = psPolynomial2DAlloc (PS_POLYNOMIAL_CHEB, nXtrend, nYtrend); - break; - - case PM_TREND_MAP: { - // binning defines the map scale relationship - psImageBinning *binning = psImageBinningAlloc(); - binning->nXruff = nXtrend; - binning->nYruff = nYtrend; - binning->nXfine = image->numCols; - binning->nYfine = image->numRows; - - trend->map = psImageMapAlloc (image, binning, stats); - psFree (binning); - break; - } - - default: - psAbort ("error"); - } - return (trend); -} - -pmTrend2D *pmTrend2DFieldAlloc (pmTrend2DMode mode, int nXfield, int nYfield, int nXtrend, int nYtrend, psStats *stats) -{ + pmTrend2D *trend = (pmTrend2D *) psAlloc(sizeof(pmTrend2D)); psMemSetDeallocator(trend, (psFreeFunc) pmTrend2DFree); @@ -102,4 +63,93 @@ // binning defines the map scale relationship psImageBinning *binning = psImageBinningAlloc(); + binning->nXruff = nXtrend; + binning->nYruff = nYtrend; + binning->nXfine = image->numCols; + binning->nYfine = image->numRows; + + trend->map = psImageMapAlloc (image, binning, stats); + psFree (binning); + break; + } + + default: + psAbort ("error"); + } + return (trend); +} + +pmTrend2D *pmTrend2DNoImageAlloc (pmTrend2DMode mode, psImageBinning *binning, psStats *stats) +{ + pmTrend2D *trend = (pmTrend2D *) psAlloc(sizeof(pmTrend2D)); + psMemSetDeallocator(trend, (psFreeFunc) pmTrend2DFree); + + trend->map = NULL; + trend->poly = NULL; + trend->stats = psMemIncrRefCounter (stats); + trend->mode = mode; + + switch (mode) { + case PM_TREND_POLY_ORD: + trend->poly = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, binning->nXruff, binning->nYruff); + // set masking somehow + for (int nx = 0; nx < trend->poly->nX + 1; nx++) { + for (int ny = 0; ny < trend->poly->nY + 1; ny++) { + if (nx + ny >= PS_MAX (trend->poly->nX, trend->poly->nY) + 1) { + trend->poly->mask[nx][ny] = 1; + } else { + trend->poly->mask[nx][ny] = 0; + } + } + } + break; + + case PM_TREND_POLY_CHEB: + trend->poly = psPolynomial2DAlloc (PS_POLYNOMIAL_CHEB, binning->nXruff, binning->nYruff); + break; + + case PM_TREND_MAP: { + // binning defines the map scale relationship + trend->map = psImageMapNoImageAlloc (binning, stats); + break; + } + + default: + psAbort ("error"); + } + return (trend); +} + +pmTrend2D *pmTrend2DFieldAlloc (pmTrend2DMode mode, int nXfield, int nYfield, int nXtrend, int nYtrend, psStats *stats) +{ + pmTrend2D *trend = (pmTrend2D *) psAlloc(sizeof(pmTrend2D)); + psMemSetDeallocator(trend, (psFreeFunc) pmTrend2DFree); + + trend->map = NULL; + trend->poly = NULL; + trend->stats = psMemIncrRefCounter (stats); + trend->mode = mode; + + switch (mode) { + case PM_TREND_POLY_ORD: + trend->poly = psPolynomial2DAlloc (PS_POLYNOMIAL_ORD, nXtrend, nYtrend); + // set masking somehow + for (int nx = 0; nx < trend->poly->nX + 1; nx++) { + for (int ny = 0; ny < trend->poly->nY + 1; ny++) { + if (nx + ny >= PS_MAX (trend->poly->nX, trend->poly->nY) + 1) { + trend->poly->mask[nx][ny] = 1; + } else { + trend->poly->mask[nx][ny] = 0; + } + } + } + break; + + case PM_TREND_POLY_CHEB: + trend->poly = psPolynomial2DAlloc (PS_POLYNOMIAL_CHEB, nXtrend, nYtrend); + break; + + case PM_TREND_MAP: { + // binning defines the map scale relationship + psImageBinning *binning = psImageBinningAlloc(); binning->nXfine = nXfield; binning->nYfine = nYfield; @@ -122,4 +172,9 @@ bool status; + assert (trend); + assert (x); + assert (y); + assert (f); + switch (trend->mode) { case PM_TREND_POLY_ORD: @@ -147,5 +202,5 @@ double result; - if (!trend) return 0.0; + assert (trend); switch (trend->mode) { @@ -168,4 +223,8 @@ psVector *result; + + assert (trend); + assert (x); + assert (y); switch (trend->mode) { Index: /branches/eam_branch_20070921/psModules/src/objects/pmTrend2D.h =================================================================== --- /branches/eam_branch_20070921/psModules/src/objects/pmTrend2D.h (revision 14988) +++ /branches/eam_branch_20070921/psModules/src/objects/pmTrend2D.h (revision 14989) @@ -5,6 +5,6 @@ * @author EAM, IfA * - * @version $Revision: 1.2.2.2 $ $Name: not supported by cvs2svn $ - * @date $Date: 2007-09-22 02:20:11 $ + * @version $Revision: 1.2.2.3 $ $Name: not supported by cvs2svn $ + * @date $Date: 2007-09-24 03:00:44 $ * Copyright 2004 Maui High Performance Computing Center, University of Hawaii */ @@ -34,4 +34,6 @@ pmTrend2D *pmTrend2DAlloc (pmTrend2DMode mode, psImage *image, int nXtrend, int nYtrend, psStats *stats); +pmTrend2D *pmTrend2DNoImageAlloc (pmTrend2DMode mode, psImageBinning *binning, psStats *stats); + // allocate a pmTrend2D tied to an abstract field with size nXfield,nYfield pmTrend2D *pmTrend2DFieldAlloc (pmTrend2DMode mode, int nXfield, int nYfield, int nXtrend, int nYtrend, psStats *stats);