IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links

Changeset 3427


Ignore:
Timestamp:
Mar 15, 2005, 1:19:36 PM (21 years ago)
Author:
Paul Price
Message:

untabified

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/doc/modules/ModulesSDRS.tex

    r3426 r3427  
    1 %%% $Id: ModulesSDRS.tex,v 1.34 2005-03-15 23:18:29 price Exp $
     1%%% $Id: ModulesSDRS.tex,v 1.35 2005-03-15 23:19:36 price Exp $
    22\documentclass[panstarrs]{panstarrs}
    33
     
    15551555\begin{verbatim}
    15561556psImage *pmCombineImages(psImage *combined, // Combined image
    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                         );
     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                        );
    15671567\end{verbatim}
    15681568
     
    16011601\begin{verbatim}
    16021602psArray *pmRejectPixels(const psArray *images, // Array of input images
    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                         );
     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                        );
    16091609\end{verbatim}
    16101610
     
    16691669
    16701670    for (int i = 0; i < nImages; i++) {
    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);
     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);
    16851685    }
    16861686
    1687     psArray *rejected = NULL;           // Array of rejected pixel lists
     1687    psArray *rejected = NULL;           // Array of rejected pixel lists
    16881688    psStats *combineStats = psStatsAlloc(PS_STAT_SAMPLE_MEAN); // Statistic to use in doing the combination
    16891689    psImage *combined = pmCombineImages(NULL, &rejected, transformed, transformedErr, transformedMask, 0,
    1690                                         NULL, numIter, sigmaClip, combineStats); // Combined image
     1690                                        NULL, numIter, sigmaClip, combineStats); // Combined image
    16911691    psArray *bad = pmRejectPixels(inputs, rejected, NULL, skyToDetector, rejThreshold, gradLimit); // Bad pix
    1692     psPixels *combinePixels = NULL;     // Pixels to combine
     1692    psPixels *combinePixels = NULL;     // Pixels to combine
    16931693    for (int i = 0; i < nImages; i++) {
    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);
     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);
    17201720    }
    17211721    psFree(bad);
     
    17231723    // Combine with no rejection
    17241724    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);
    17261726    psFree(combineStats);
    17271727    psFree(combinePixels);
     
    17631763\begin{verbatim}
    17641764typedef enum {
    1765     PM_SUBTRACTION_KERNEL_POIS,         // POIS kernel --- delta functions
    1766     PM_SUBTRACTION_KERNEL_ISIS          // ISIS kernel --- gaussians modified by polynomials
     1765    PM_SUBTRACTION_KERNEL_POIS,         // POIS kernel --- delta functions
     1766    PM_SUBTRACTION_KERNEL_ISIS          // ISIS kernel --- gaussians modified by polynomials
    17671767} pmSubtractionKernelsType;
    17681768\end{verbatim}
     
    17741774\begin{verbatim}
    17751775typedef struct {
    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)
     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)
    17801780    int subIndex;                       // Index of kernel to be subtracted to maintain flux conservation
    1781     psArray *preCalc;                   // Array of images containing pre-calculated kernel (to
     1781    psArray *preCalc;                   // Array of images containing pre-calculated kernel (to
    17821782                                        // accelerate ISIS; don't use for POIS)
    17831783} pmSubtractionKernels;
     
    18461846\begin{verbatim}
    18471847typedef struct {
    1848     int x, y;                           // Position
    1849     psImage *matrix;                    // Associated matrix
    1850     psVector *vector;                   // Associated vector
    1851     psStampStatus status;               // Status of stamp
     1848    int x, y;                           // Position
     1849    psImage *matrix;                    // Associated matrix
     1850    psVector *vector;                   // Associated vector
     1851    psStampStatus status;               // Status of stamp
    18521852} pmStamp;
    18531853\end{verbatim}
     
    18671867\begin{verbatim}
    18681868typedef enum {
    1869     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
     1869    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
    18731873} pmStampStatus;
    18741874\end{verbatim}
     
    18761876\begin{verbatim}
    18771877psArary *pmSubtractionFindStamps(psArray *stamps, // Output stamps, or NULL
    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                                 );
     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                                );
    18851885\end{verbatim}
    18861886
     
    19131913\begin{verbatim}
    19141914bool pmSubtractionCalculateEquation(psArray *stamps, // The stamps for which to calculate the equation
    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                                     );
     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                                    );
    19201920\end{verbatim}
    19211921
     
    19591959
    19601960\begin{verbatim}
    1961 psVector *pmSubtractionSolveEquation(psVector *solution,        // Solution vector, or NULL
    1962                                      const psArray *stamps // Array of stamps
    1963                                      );
     1961psVector *pmSubtractionSolveEquation(psVector *solution,        // Solution vector, or NULL
     1962                                     const psArray *stamps // Array of stamps
     1963                                     );
    19641964\end{verbatim}
    19651965
     
    19811981\begin{verbatim}
    19821982bool pmSubtractionRejectStamps(psArray *stamps, // Array of stamps to check for rejection
    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                                );
     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                               );
    19921992\end{verbatim}
    19931993
     
    20592059
    20602060    // Generate kernel basis functions
    2061     psArray *kernels = NULL;            // Array of kernel basis functions
     2061    psArray *kernels = NULL;            // Array of kernel basis functions
    20622062    switch (kernelType) {
    20632063      case PM_SUBTRACTION_KERNEL_POIS:
    2064         // Create the kernel basis functions
    2065         kernels = pmSubtractionKernelsGeneratePOIS(kernelHalfSize, spatialOrder);
    2066         break;
     2064        // Create the kernel basis functions
     2065        kernels = pmSubtractionKernelsGeneratePOIS(kernelHalfSize, spatialOrder);
     2066        break;
    20672067      case PM_SUBTRACTION_KERNEL_ISIS:
    2068         kernels = pmSubtractionKernelsGenerateISIS(sigmas, polyOrders, kernelHalfSize, spatialOrder);
    2069         break;
     2068        kernels = pmSubtractionKernelsGenerateISIS(sigmas, polyOrders, kernelHalfSize, spatialOrder);
     2069        break;
    20702070      default:
    2071         barf();
     2071        barf();
    20722072    }
    20732073
    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?
     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?
    20772077
    20782078    // Iterate for a solution
    20792079    for (int iter = 0; iter < numIter && rejected; iter++) {
    20802080
    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);
     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);
    20922092    }
    20932093
     
    21402140
    21412141  keep = psRegionAlloc (XBORDER, image->nCol - YBORDER,
    2142                         YBORDER, image->nRow - YBORDER);
     2142                        YBORDER, image->nRow - YBORDER);
    21432143 
    21442144  Sky = stats->median;
Note: See TracChangeset for help on using the changeset viewer.