IPP Software Navigation Tools IPP Links Communication Pan-STARRS Links
wiki:Stack_Algorithm

Version 1 (modified by watersc1, 13 years ago) ( diff )

--

Intro

Input selection

PSF Matching

Pixel Rejection

The following conventions in the remainder of this discussion. Lower case values are for a single pixel on image "i", upper case values are values that are constant for a given image "i" (i.e., all pixels on image i have that parameter), greek values are combined values (the stacked values).

# Set up For a given pixel, define four parameters: f_i pixel value v_i pixel variance value m_i input mask value, 0 = bad, 1 = good s_i suspect mask value, 0 = not suspect, 1 = suspect W_i image weight

There are also a set of configurable options: NMIN = 5 Minimum number of values for KMM test MIN_UNI = 0.05 Minimum probability to accept a bimodal solution SYS_FRAC = 0.1 Systematic variance fraction NSIG_REJ = 3.5 N-sigma rejection threshold DISCARD = 0.2 Fraction of values to exclude in Olympic mean ITER = 0.5 Number of rejection iterations per input

After convolution to the target PSF, we add an additional variance (ADDV_i) term equal to the chi2 value found during the PSF matching. For 1d Gaussian convolution, we have no real concept of this value, so all chi2 (and ADDV_i) are defined to be 1.0;

v_i = v_i + ADDV_i

The image weights are calculated by taking the inverse of the robust median of the variance image:

W_i ~ 1 / sum_{x,y} v_i_{x,y}

# Calulate initial stack values

  • Calculate the weighted mean:

mu = sum_i( f_i * W_i) / sum_i(W_i) sigma = 1 / sum_i(1 / W_i) nu = sum_i(m_i)

# Rejection loop

This loop is repeated ITER * Ninput times (minimum of once).

  • Determine distribution values:

if (nu > NMIN) do Mixture Model test:

mu_KMM mean of most popular gaussian mode sigma_KMM sigma of most popular gaussian mode pi_KMM fraction of inputs in most popular gaussian mode rho_UNI Probability that inputs are unimodal

  • Set rejection limits for each input:

if (nu > NMIN) sys_var = sigma_KMM2 else sys_var = SYS_FRAC * f_i

limit_i = (NSIG_REJ)2 * (v_i + sys_var)

  • Determine distribution median

if (nu > NMIN) mu_median = mu_KMM else mu_median = olympic_mean where: olympic_mean = sum_i(f_i * W_i) / sum_i(W_i) from int((DISCARD * nu + 0.5)/2.0) to int((DISCARD * nu + 0.5)/2.0) + int(nu) - int(DISCARD * nu + 0.5)

  • Determine most discrepant point

if ( (f_i - mu_median)2 > limit_i ) AND

( (f_i - mu_median)2 / limit_i > delta ) THEN

delta = (f_i - mu_median)2 / limit_i delta_index = i

  • Do rejection. If we have suspect pixels, throw those out first (so they'll be skipped in the next iteration). Otherwise, remove the most discrepant point

if DEFINED(dev_i)

if (s_i != 0) m_i = 0 else if (i = delta_index) m_i = 0

  • End rejection loop.

This completes the initial stack stage of ppStack. We then use the set of rejected pixels to identify the source pixels in the unconvolved input images. These are masked out, and the images are reconvolved and a final stack is constructed using the same algorithm, however the rejection loop is skipped (so the weighted mean is all that is performed).

Output Products

Note: See TracWiki for help on using the wiki.