Index: trunk/psModules/src/imcombine/pmSubtractionKernels.h =================================================================== --- trunk/psModules/src/imcombine/pmSubtractionKernels.h (revision 26157) +++ trunk/psModules/src/imcombine/pmSubtractionKernels.h (revision 26893) @@ -10,4 +10,7 @@ PM_SUBTRACTION_KERNEL_POIS, ///< Pan-STARRS Optimal Image Subtraction --- delta functions PM_SUBTRACTION_KERNEL_ISIS, ///< Traditional kernel --- gaussians modified by polynomials + PM_SUBTRACTION_KERNEL_ISIS_RADIAL, ///< ISIS + higher-order radial Hermitians + PM_SUBTRACTION_KERNEL_HERM, ///< Hermitian polynomial kernels + PM_SUBTRACTION_KERNEL_DECONV_HERM, ///< Deconvolved Hermitian polynomial kernels PM_SUBTRACTION_KERNEL_SPAM, ///< Summed Pixels for Advanced Matching --- summed delta functions PM_SUBTRACTION_KERNEL_FRIES, ///< Fibonacci Radius Increases Excellence of Subtraction @@ -29,9 +32,10 @@ pmSubtractionKernelsType type; ///< Type of kernels --- allowing the use of multiple kernels psString description; ///< Description of the kernel parameters + int xMin, xMax, yMin, yMax; ///< Bounds of image (for normalisation) long num; ///< Number of kernel components (not including the spatial ones) - psVector *u, *v; ///< Offset (for POIS) or polynomial order (for ISIS) - psVector *widths; ///< Gaussian FWHMs (ISIS) + psVector *u, *v; ///< Offset (for POIS) or polynomial order (for ISIS, HERM or DECONV_HERM) + psVector *widths; ///< Gaussian FWHMs (ISIS, HERM or DECONV_HERM) psVector *uStop, *vStop; ///< Width of kernel element (SPAM,FRIES only) - psArray *preCalc; ///< Array of images containing pre-calculated kernel (for ISIS) + psArray *preCalc; ///< Array of images containing pre-calculated kernel (for ISIS, HERM or DECONV_HERM) float penalty; ///< Penalty for wideness psVector *penalties; ///< Penalty for each kernel component @@ -41,10 +45,27 @@ int bgOrder; ///< The order for the background fitting pmSubtractionMode mode; ///< Mode for subtraction - int numCols, numRows; ///< Size of image (for normalisation), or zero to use image provided psVector *solution1, *solution2; ///< Solution for the PSF matching // Quality information float mean, rms; ///< Mean and RMS of chi^2 from stamps int numStamps; ///< Number of good stamps + float fSigResMean; ///< mean fractional stdev of residuals + float fSigResStdev; ///< stdev of fractional stdev of residuals + float fMaxResMean; ///< mean fractional positive swing in residuals + float fMaxResStdev; ///< stdev of fractional positive swing in residuals + float fMinResMean; ///< mean fractional negative swing in residuals + float fMinResStdev; ///< stdev of fractional negative swing in residuals + psArray *sampleStamps; ///< array of brightest set of stamps for output visualizations } pmSubtractionKernels; + +// pmSubtractionKernels->preCalc is an array of pmSubtractionKernelPreCalc structures +typedef struct { + psVector *uCoords; // used by RINGS + psVector *vCoords; // used by RINGS + psVector *poly; // used by RINGS + + psVector *xKernel; // used by ISIS, HERM, DECONV_HERM + psVector *yKernel; // used by ISIS, HERM, DECONV_HERM + psKernel *kernel; // used by ISIS, HERM, DECONV_HERM +} pmSubtractionKernelPreCalc; // Assertion to check pmSubtractionKernels @@ -64,4 +85,19 @@ PS_ASSERT_VECTOR_SIZE((KERNELS)->widths, (KERNELS)->num, RETURNVALUE); \ } \ + if ((KERNELS)->type == PM_SUBTRACTION_KERNEL_ISIS_RADIAL) { \ + PS_ASSERT_VECTOR_NON_NULL((KERNELS)->widths, RETURNVALUE); \ + PS_ASSERT_VECTOR_TYPE((KERNELS)->widths, PS_TYPE_F32, RETURNVALUE); \ + PS_ASSERT_VECTOR_SIZE((KERNELS)->widths, (KERNELS)->num, RETURNVALUE); \ + } \ + if ((KERNELS)->type == PM_SUBTRACTION_KERNEL_HERM) { \ + PS_ASSERT_VECTOR_NON_NULL((KERNELS)->widths, RETURNVALUE); \ + PS_ASSERT_VECTOR_TYPE((KERNELS)->widths, PS_TYPE_F32, RETURNVALUE); \ + PS_ASSERT_VECTOR_SIZE((KERNELS)->widths, (KERNELS)->num, RETURNVALUE); \ + } \ + if ((KERNELS)->type == PM_SUBTRACTION_KERNEL_DECONV_HERM) { \ + PS_ASSERT_VECTOR_NON_NULL((KERNELS)->widths, RETURNVALUE); \ + PS_ASSERT_VECTOR_TYPE((KERNELS)->widths, PS_TYPE_F32, RETURNVALUE); \ + PS_ASSERT_VECTOR_SIZE((KERNELS)->widths, (KERNELS)->num, RETURNVALUE); \ + } \ if ((KERNELS)->uStop || (KERNELS)->vStop) { \ PS_ASSERT_VECTOR_NON_NULL((KERNELS)->uStop, RETURNVALUE); \ @@ -99,4 +135,16 @@ } +// Generate 1D convolution kernel for ISIS +psVector *pmSubtractionKernelISIS(float sigma, // Gaussian width + int order, // Polynomial order + int size // Kernel half-size + ); + +// Generate 1D convolution kernel for HERM (normalized for 2D) +psVector *pmSubtractionKernelHERM(float sigma, // Gaussian width + int order, // Polynomial order + int size // Kernel half-size + ); + /// Generate a delta-function grid for subtraction kernels (like the POIS kernel) bool p_pmSubtractionKernelsAddGrid(pmSubtractionKernels *kernels, ///< The subtraction kernels to append to @@ -114,6 +162,17 @@ int spatialOrder, ///< Order of spatial variations float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity pmSubtractionMode mode ///< Mode for subtraction ); + +/// Allocator for pre-calculated kernel data structure +pmSubtractionKernelPreCalc *pmSubtractionKernelPreCalcAlloc( + pmSubtractionKernelsType type, ///< type of kernel to allocate (not all can be pre-calculated) + int uOrder, ///< order in x-direction + int vOrder, ///< order in x-direction + int size, ///< Half-size of the kernel + float sigma ///< sigma of gaussian kernel + ); + /// Generate POIS kernels @@ -121,4 +180,5 @@ int spatialOrder, ///< Order of spatial variations float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity pmSubtractionMode mode ///< Mode for subtraction ); @@ -130,4 +190,5 @@ const psVector *orders, ///< Polynomial order of gaussians float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity pmSubtractionMode mode ///< Mode for subtraction ); @@ -139,6 +200,37 @@ const psVector *orders, ///< Polynomial order of gaussians float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity pmSubtractionMode mode ///< Mode for subtraction ); + +/// Generate ISIS + RADIAL_HERM kernels +pmSubtractionKernels *pmSubtractionKernelsISIS_RADIAL(int size, ///< Half-size of the kernel + int spatialOrder, ///< Order of spatial variations + const psVector *fwhms, ///< Gaussian FWHMs + const psVector *orders, ///< Polynomial order of gaussians + float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity + pmSubtractionMode mode ///< Mode for subtraction + ); + +/// Generate HERM kernels +pmSubtractionKernels *pmSubtractionKernelsHERM(int size, ///< Half-size of the kernel + int spatialOrder, ///< Order of spatial variations + const psVector *fwhms, ///< Gaussian FWHMs + const psVector *orders, ///< order of hermitian polynomials + float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity + pmSubtractionMode mode ///< Mode for subtraction + ); + +/// Generate DECONV_HERM kernels +pmSubtractionKernels *pmSubtractionKernelsDECONV_HERM(int size, ///< Half-size of the kernel + int spatialOrder, ///< Order of spatial variations + const psVector *fwhms, ///< Gaussian FWHMs + const psVector *orders, ///< order of hermitian polynomials + float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity + pmSubtractionMode mode ///< Mode for subtraction + ); /// Generate SPAM kernels @@ -148,4 +240,5 @@ int binning, ///< Kernel binning factor float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity pmSubtractionMode mode ///< Mode for subtraction ); @@ -156,4 +249,5 @@ int inner, ///< Inner radius to preserve unbinned float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity pmSubtractionMode mode ///< Mode for subtraction ); @@ -166,4 +260,5 @@ int inner, ///< Inner radius containing grid of delta functions float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity pmSubtractionMode mode ///< Mode for subtraction ); @@ -175,4 +270,5 @@ int ringsOrder, ///< Polynomial order float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity pmSubtractionMode mode ///< Mode for subtraction ); @@ -189,4 +285,5 @@ int ringsOrder, ///< Polynomial order for RINGS float penalty, ///< Penalty for wideness + psRegion bounds, ///< Bounds for validity pmSubtractionMode mode ///< Mode for subtraction ); @@ -196,4 +293,5 @@ const char *description, ///< Description of kernel int bgOrder, ///< Polynomial order for background fitting + psRegion bounds, ///< Bounds for validity pmSubtractionMode mode ///< Mode for subtraction );