Changeset 3427
- Timestamp:
- Mar 15, 2005, 1:19:36 PM (21 years ago)
- File:
-
- 1 edited
-
trunk/doc/modules/ModulesSDRS.tex (modified) (15 diffs)
Legend:
- Unmodified
- Added
- Removed
-
trunk/doc/modules/ModulesSDRS.tex
r3426 r3427 1 %%% $Id: ModulesSDRS.tex,v 1.3 4 2005-03-15 23:18:29price Exp $1 %%% $Id: ModulesSDRS.tex,v 1.35 2005-03-15 23:19:36 price Exp $ 2 2 \documentclass[panstarrs]{panstarrs} 3 3 … … 1555 1555 \begin{verbatim} 1556 1556 psImage *pmCombineImages(psImage *combined, // Combined image 1557 psArray **questionablePixels, // Array of rejection masks1558 const psArray *images, // Array of input images1559 const psArray *errors, // Array of input error images1560 const psArray *masks,// Array of input masks1561 unsigned int maskVal, // Mask value1562 const psArray *pixels,// Pixels to combine1563 int numIter,// Number of rejection iterations1564 float sigmaClip, // Number of standard deviations at which to reject1565 const psStats *stats // Statistics to use in the combination1566 );1557 psArray **questionablePixels, // Array of rejection masks 1558 const psArray *images, // Array of input images 1559 const psArray *errors, // Array of input error images 1560 const psArray *masks,// Array of input masks 1561 unsigned int maskVal, // Mask value 1562 const psArray *pixels, // Pixels to combine 1563 int numIter, // Number of rejection iterations 1564 float sigmaClip, // Number of standard deviations at which to reject 1565 const psStats *stats // Statistics to use in the combination 1566 ); 1567 1567 \end{verbatim} 1568 1568 … … 1601 1601 \begin{verbatim} 1602 1602 psArray *pmRejectPixels(const psArray *images, // Array of input images 1603 const psArray *pixels, // These are the pixels which were rejected in the combination1604 const psArray *inToOut, // Transformations from input to output system1605 const psArray *outToIn, // Transformations from output to input system1606 float rejThreshold, // Rejection threshold1607 float gradLimit // Gradient limit1608 );1603 const psArray *pixels, // These are the pixels which were rejected in the combination 1604 const psArray *inToOut, // Transformations from input to output system 1605 const psArray *outToIn, // Transformations from output to input system 1606 float rejThreshold, // Rejection threshold 1607 float gradLimit // Gradient limit 1608 ); 1609 1609 \end{verbatim} 1610 1610 … … 1669 1669 1670 1670 for (int i = 0; i < nImages; i++) { 1671 psPixels *blanks = NULL;// List of blank pixels1672 transformed->data[i] = psImageTransform(NULL, &blanks, inputs->data[i],1673 inputMask->data[i], inputMaskVal, NAN, skyToDetector,1674 combineRegion, NULL, PS_INTERPOLATE_BILINEAR);1675 transformedErr->data[i] = psImageTransform(NULL, NULL, inputsErr->data[i], inputMask->data[i],1676 inputMaskVal, NAN, skyToDetector, combineRegion, NULL,1677 PS_INTERPOLATE_BILINEAR_VARIANCE);1678 psImage *skyImage = transformed->data[i]; // Dereference the transformed image1679 psRegion *blankRegion = psRegionAlloc(0, 0, skyImage->numCols, skyImage->numRows); // Size of1680 // transformed1681 // image1682 transformedMask->data[i] = psPixelsToMask(NULL, blanks, blankRegion, PS_MASK_BLANK);1683 psFree(blankRegion);1684 psFree(blanks);1671 psPixels *blanks = NULL; // List of blank pixels 1672 transformed->data[i] = psImageTransform(NULL, &blanks, inputs->data[i], 1673 inputMask->data[i], inputMaskVal, NAN, skyToDetector, 1674 combineRegion, NULL, PS_INTERPOLATE_BILINEAR); 1675 transformedErr->data[i] = psImageTransform(NULL, NULL, inputsErr->data[i], inputMask->data[i], 1676 inputMaskVal, NAN, skyToDetector, combineRegion, NULL, 1677 PS_INTERPOLATE_BILINEAR_VARIANCE); 1678 psImage *skyImage = transformed->data[i]; // Dereference the transformed image 1679 psRegion *blankRegion = psRegionAlloc(0, 0, skyImage->numCols, skyImage->numRows); // Size of 1680 // transformed 1681 // image 1682 transformedMask->data[i] = psPixelsToMask(NULL, blanks, blankRegion, PS_MASK_BLANK); 1683 psFree(blankRegion); 1684 psFree(blanks); 1685 1685 } 1686 1686 1687 psArray *rejected = NULL; // Array of rejected pixel lists1687 psArray *rejected = NULL; // Array of rejected pixel lists 1688 1688 psStats *combineStats = psStatsAlloc(PS_STAT_SAMPLE_MEAN); // Statistic to use in doing the combination 1689 1689 psImage *combined = pmCombineImages(NULL, &rejected, transformed, transformedErr, transformedMask, 0, 1690 NULL, numIter, sigmaClip, combineStats); // Combined image1690 NULL, numIter, sigmaClip, combineStats); // Combined image 1691 1691 psArray *bad = pmRejectPixels(inputs, rejected, NULL, skyToDetector, rejThreshold, gradLimit); // Bad pix 1692 psPixels *combinePixels = NULL; // Pixels to combine1692 psPixels *combinePixels = NULL; // Pixels to combine 1693 1693 for (int i = 0; i < nImages; i++) { 1694 psPixels *badSource = psPixelsTransform(NULL, bad->data[i], skyToDetector); // Bad pixels on the input1695 psImage *badMask = psPixelsToMask(NULL, badSource, PS_MASK_COSMICRAY); // Mask image for the input1696 (void)psBinaryOp(inputMask->data[i], inputMask->data[i], "|", badMask);// Put CRs into original mask1697 psFree(badSource);1698 psFree(badMask);1699 1700 combinePixels = psPixelsConcatenate(redo, bad->data[i]);1701 1702 // Update transformed image1703 psPixels *blanks = NULL;// List of blank pixels1704 transformed->data[i] = psImageTransform(transformed->data[i], &blanks, inputs->data[i],1705 inputMask->data[i], inputMaskVal | PS_MASK_COSMICRAY, NAN,1706 skyToDetector, combineRegion, bad->data[i],1707 PS_INTERPOLATE_BILINEAR);1708 transformedErr->data[i] = psImageTransform(transformedErr->data[i], NULL, inputsErr->data[i],1709 inputMask->data[i], inputMaskVal | PS_MASK_COSMICRAY,1710 NAN, skyToDetector, combineRegion, bad->data[i],1711 PS_INTERPOLATE_BILINEAR_VARIANCE);1712 psImage *skyImage = transformed->data[i]; // Dereference the transformed image1713 psRegion *blankRegion = psRegionAlloc(0, 0, skyImage->numCols, skyImage->numRows); // Size of1714 // transformed1715 // image1716 transformedMask->data[i] = psPixelsToMask(transformedMask->data[i], blanks, blankRegion,1717 PS_MASK_BLANK);1718 psFree(blankRegion);1719 psFree(blanks);1694 psPixels *badSource = psPixelsTransform(NULL, bad->data[i], skyToDetector); // Bad pixels on the input 1695 psImage *badMask = psPixelsToMask(NULL, badSource, PS_MASK_COSMICRAY); // Mask image for the input 1696 (void)psBinaryOp(inputMask->data[i], inputMask->data[i], "|", badMask); // Put CRs into original mask 1697 psFree(badSource); 1698 psFree(badMask); 1699 1700 combinePixels = psPixelsConcatenate(redo, bad->data[i]); 1701 1702 // Update transformed image 1703 psPixels *blanks = NULL; // List of blank pixels 1704 transformed->data[i] = psImageTransform(transformed->data[i], &blanks, inputs->data[i], 1705 inputMask->data[i], inputMaskVal | PS_MASK_COSMICRAY, NAN, 1706 skyToDetector, combineRegion, bad->data[i], 1707 PS_INTERPOLATE_BILINEAR); 1708 transformedErr->data[i] = psImageTransform(transformedErr->data[i], NULL, inputsErr->data[i], 1709 inputMask->data[i], inputMaskVal | PS_MASK_COSMICRAY, 1710 NAN, skyToDetector, combineRegion, bad->data[i], 1711 PS_INTERPOLATE_BILINEAR_VARIANCE); 1712 psImage *skyImage = transformed->data[i]; // Dereference the transformed image 1713 psRegion *blankRegion = psRegionAlloc(0, 0, skyImage->numCols, skyImage->numRows); // Size of 1714 // transformed 1715 // image 1716 transformedMask->data[i] = psPixelsToMask(transformedMask->data[i], blanks, blankRegion, 1717 PS_MASK_BLANK); 1718 psFree(blankRegion); 1719 psFree(blanks); 1720 1720 } 1721 1721 psFree(bad); … … 1723 1723 // Combine with no rejection 1724 1724 combined = pmCombineImages(combined, NULL, transformed, transformedErr, transformedMask, 1725 PS_MASK_BLANK, combinePixels, 0, 0.0, combineStats);1725 PS_MASK_BLANK, combinePixels, 0, 0.0, combineStats); 1726 1726 psFree(combineStats); 1727 1727 psFree(combinePixels); … … 1763 1763 \begin{verbatim} 1764 1764 typedef enum { 1765 PM_SUBTRACTION_KERNEL_POIS, // POIS kernel --- delta functions1766 PM_SUBTRACTION_KERNEL_ISIS // ISIS kernel --- gaussians modified by polynomials1765 PM_SUBTRACTION_KERNEL_POIS, // POIS kernel --- delta functions 1766 PM_SUBTRACTION_KERNEL_ISIS // ISIS kernel --- gaussians modified by polynomials 1767 1767 } pmSubtractionKernelsType; 1768 1768 \end{verbatim} … … 1774 1774 \begin{verbatim} 1775 1775 typedef struct { 1776 pmSubtractionKernelType type; // Type of kernels --- allowing the use of multiple kernels1777 psVector *u, *v; // Offset (for POIS) or polynomial order (for ISIS)1778 psVector *sigma; // Width of Gaussian (for ISIS)1779 psVector *xOrder, *yOrder; // Spatial polynomial order (for all)1776 pmSubtractionKernelType type; // Type of kernels --- allowing the use of multiple kernels 1777 psVector *u, *v; // Offset (for POIS) or polynomial order (for ISIS) 1778 psVector *sigma; // Width of Gaussian (for ISIS) 1779 psVector *xOrder, *yOrder; // Spatial polynomial order (for all) 1780 1780 int subIndex; // Index of kernel to be subtracted to maintain flux conservation 1781 psArray *preCalc; // Array of images containing pre-calculated kernel (to1781 psArray *preCalc; // Array of images containing pre-calculated kernel (to 1782 1782 // accelerate ISIS; don't use for POIS) 1783 1783 } pmSubtractionKernels; … … 1846 1846 \begin{verbatim} 1847 1847 typedef struct { 1848 int x, y; // Position1849 psImage *matrix; // Associated matrix1850 psVector *vector; // Associated vector1851 psStampStatus status; // Status of stamp1848 int x, y; // Position 1849 psImage *matrix; // Associated matrix 1850 psVector *vector; // Associated vector 1851 psStampStatus status; // Status of stamp 1852 1852 } pmStamp; 1853 1853 \end{verbatim} … … 1867 1867 \begin{verbatim} 1868 1868 typedef enum { 1869 PM_STAMP_USED, // Use this stamp1870 PM_STAMP_REJECTED, // This stamp has been rejected1871 PM_STAMP_RECALC, // Having been reset, this stamp needs to be recalculated1872 PM_STAMP_NONE // No stamp in this region1869 PM_STAMP_USED, // Use this stamp 1870 PM_STAMP_REJECTED, // This stamp has been rejected 1871 PM_STAMP_RECALC, // Having been reset, this stamp needs to be recalculated 1872 PM_STAMP_NONE // No stamp in this region 1873 1873 } pmStampStatus; 1874 1874 \end{verbatim} … … 1876 1876 \begin{verbatim} 1877 1877 psArary *pmSubtractionFindStamps(psArray *stamps, // Output stamps, or NULL 1878 const psImage *image, // Image for which to find stamps1879 const psImage *mask, // Mask1880 unsigned int maskVal, // Value for mask1881 float threshold, // Threshold for stamps in the image1882 int xNum, int yNum, // Number of stamps in x and y1883 int border // Border around image to ignore (should be size of kernel)1884 );1878 const psImage *image, // Image for which to find stamps 1879 const psImage *mask, // Mask 1880 unsigned int maskVal, // Value for mask 1881 float threshold, // Threshold for stamps in the image 1882 int xNum, int yNum, // Number of stamps in x and y 1883 int border // Border around image to ignore (should be size of kernel) 1884 ); 1885 1885 \end{verbatim} 1886 1886 … … 1913 1913 \begin{verbatim} 1914 1914 bool pmSubtractionCalculateEquation(psArray *stamps, // The stamps for which to calculate the equation 1915 const psImage *reference, // Reference image1916 const psImage *input, // Input image1917 const psSubtractionKernels *kernels, // The kernel basis functions1918 int footprint // Half-size of region over which to calculate equation1919 );1915 const psImage *reference, // Reference image 1916 const psImage *input, // Input image 1917 const psSubtractionKernels *kernels, // The kernel basis functions 1918 int footprint // Half-size of region over which to calculate equation 1919 ); 1920 1920 \end{verbatim} 1921 1921 … … 1959 1959 1960 1960 \begin{verbatim} 1961 psVector *pmSubtractionSolveEquation(psVector *solution, // Solution vector, or NULL1962 const psArray *stamps // Array of stamps1963 );1961 psVector *pmSubtractionSolveEquation(psVector *solution, // Solution vector, or NULL 1962 const psArray *stamps // Array of stamps 1963 ); 1964 1964 \end{verbatim} 1965 1965 … … 1981 1981 \begin{verbatim} 1982 1982 bool pmSubtractionRejectStamps(psArray *stamps, // Array of stamps to check for rejection 1983 psImage *mask, // Mask image1984 unsigned int badStampMaskVal, // Value to use in mask for bad stamp1985 int footprint, // Region to mask if stamp is bad1986 float sigmaRej, // Number of RMS deviations above zero at which to reject1987 const psImage *refImage, // Reference image1988 const psImage *inImage, // Input image1989 const psVector *solution, // Solution vector1990 const pmSubtractionKernels *kernels // Array of kernel parameters1991 );1983 psImage *mask, // Mask image 1984 unsigned int badStampMaskVal, // Value to use in mask for bad stamp 1985 int footprint, // Region to mask if stamp is bad 1986 float sigmaRej, // Number of RMS deviations above zero at which to reject 1987 const psImage *refImage, // Reference image 1988 const psImage *inImage, // Input image 1989 const psVector *solution, // Solution vector 1990 const pmSubtractionKernels *kernels // Array of kernel parameters 1991 ); 1992 1992 \end{verbatim} 1993 1993 … … 2059 2059 2060 2060 // Generate kernel basis functions 2061 psArray *kernels = NULL; // Array of kernel basis functions2061 psArray *kernels = NULL; // Array of kernel basis functions 2062 2062 switch (kernelType) { 2063 2063 case PM_SUBTRACTION_KERNEL_POIS: 2064 // Create the kernel basis functions2065 kernels = pmSubtractionKernelsGeneratePOIS(kernelHalfSize, spatialOrder);2066 break;2064 // Create the kernel basis functions 2065 kernels = pmSubtractionKernelsGeneratePOIS(kernelHalfSize, spatialOrder); 2066 break; 2067 2067 case PM_SUBTRACTION_KERNEL_ISIS: 2068 kernels = pmSubtractionKernelsGenerateISIS(sigmas, polyOrders, kernelHalfSize, spatialOrder);2069 break;2068 kernels = pmSubtractionKernelsGenerateISIS(sigmas, polyOrders, kernelHalfSize, spatialOrder); 2069 break; 2070 2070 default: 2071 barf();2071 barf(); 2072 2072 } 2073 2073 2074 psArray *stamps = NULL; // Array of stamps2075 psVector *kernelCoeffs = NULL; // Coefficients for the kernels2076 bool rejected = true; // Did we reject a stamp in the last iteration?2074 psArray *stamps = NULL; // Array of stamps 2075 psVector *kernelCoeffs = NULL; // Coefficients for the kernels 2076 bool rejected = true; // Did we reject a stamp in the last iteration? 2077 2077 2078 2078 // Iterate for a solution 2079 2079 for (int iter = 0; iter < numIter && rejected; iter++) { 2080 2080 2081 // Find stamps2082 stamps = pmSubtractionFindStamps(stamps, reference, subMask, maskVal | PS_MASK_BAD_STAMP,2083 stampThreshold, nStampsX, nStampsY, stampSize, kernelHalfSize);2084 2085 // Generate and solve matrix equations2086 (void)pmSubtractionCalculateEquation(stamps, reference, input, kernels, stampSize);2087 kernelCoeffs = pmSubtractionSolveEquation(kernelCoeffs, stamps);2088 2089 // Reject bad stamps2090 rejected = pmSubtractionRejectStamps(stamps, subMask, PS_MASK_BAD_STAMP, stampSize, sigmaRej,2091 reference, input, kernelCoeffs, kernels);2081 // Find stamps 2082 stamps = pmSubtractionFindStamps(stamps, reference, subMask, maskVal | PS_MASK_BAD_STAMP, 2083 stampThreshold, nStampsX, nStampsY, stampSize, kernelHalfSize); 2084 2085 // Generate and solve matrix equations 2086 (void)pmSubtractionCalculateEquation(stamps, reference, input, kernels, stampSize); 2087 kernelCoeffs = pmSubtractionSolveEquation(kernelCoeffs, stamps); 2088 2089 // Reject bad stamps 2090 rejected = pmSubtractionRejectStamps(stamps, subMask, PS_MASK_BAD_STAMP, stampSize, sigmaRej, 2091 reference, input, kernelCoeffs, kernels); 2092 2092 } 2093 2093 … … 2140 2140 2141 2141 keep = psRegionAlloc (XBORDER, image->nCol - YBORDER, 2142 YBORDER, image->nRow - YBORDER);2142 YBORDER, image->nRow - YBORDER); 2143 2143 2144 2144 Sky = stats->median;
Note:
See TracChangeset
for help on using the changeset viewer.
