psReadoutCombine() questions

[gusciora@…]

What are the valid types for parameters:

output
inputs
zero
scale

What are the valid ranges for input parameters:

gain
readnoise

The SDR states that if the input parameters are not of the same size, then an
error should be stated. It also states that input images need not be of the
same size. Which input parameters must be of the same size? How should the
offsets be used in determining whether the input images are of the correct
size?

How should the zero vector be added to the input images? Does this vector
have one element per input image, or per row, or per column? Does this
addition happen after the pixels outside (fracLow, fracHigh) are rejected. or
before?

How should the scale vector multiply the input images? Does this vector have
one element per input image, or per row, or per column? Does this multiply
happen after the pixels outside (fracLow, fracHigh) are rejected. or before?

How are the various offsets in the input readouts, and the input images used,
along with pixel coordinates in those images, be used in determing which
!es correspond to which pixels in putout images? What is the formula?

How is the noise to be used?

[gusciora@…]

I have an algorithmic question. With the various offsets in the output image,
the input readout, and the input image, it seems that a specific pixel (i, j) in
the output image can map to an arbitrary position (x, y) in one of the input
images. In order to combine all of these input, I envision creating a psVector
for every pixel in the output image, with a length equal to the number of input
readouts, then stepping through the input readouts in order to fill those
psVectors, then applying the psVectorStats function to them, then creating the
actual output image.

Is this the basic idea that you had in mind? Or, am I misinterpreting
something? This could be very memory intensive, if there are many input images.

[gusciora@…]

I implemented this module and for the unresolved questions, I coded I thought
was the most logical choice. The basic algorithm is this:

1: Loop through all input readouts. Keep track of the max/min row/column

coordinates in output image space: the max row is

max(psReadout->row0 + psReadout->image->row0 +

psReadout->image->numRows)

similar calculation for the columns. Am I using the ->row0 member

correctly?

2: I compare max/min row/col to the output image and produce an error if

they are outside the range of that output.

3: I assume that the psVectors zero and scale must have as many elements

as the number of input psReadouts. If zero/scale have less elements, I
generate an error and return NULL. If more, I generate a warning and
proceed.

4: Generate an error if more that 1 stats options is requested.

5: I generate a single psVector with a size equal to the number of input

psReadouts. I loop through every pixel in the output image. For that
pixel, I loop through all input readouts and determine if they have a
pixel which corresponds to the output pixel. If so I save it in that
psVector. Then I save all the pixels between fracLow and fracHigh,
assuming there "nKeep" of them. If the zero vector is non NULL, I add
the i-th element of the zero vector to the pixels from the i-th
readout. I multiply the i-th element of the scale vector by the i-th
readout. If gain>0.0, I multiply that by all pixels in this vector.
If readnoise is non-negative, I add it to all pixels. Then I generate
the requested statistical quantity from that psVector. I set the
output image pixel to be the result of that statistic.

Let me know if any of the above is not correct.

[Paul Price]

Valid types for output, inputs, zero and scale are S16, S32 and F32.

If gain is positive (i.e., > 0), or the readnoise is non-negative (i.e., >= 0),
then they are used to provide Poisson weights [ variance = (gain*value +
readnoise*readnoise)/gain ].

All the lists or vectors must be of the same length, but the images (specified
within the lists or vectors) may have different sizes.

The zero and scale vectors are applied to the appropriate image before any
rejection is performed.

The offsets in the inputs are used in a similar way as the flat-fielding, except
that here there are multiple offsets --- it's a matter of working out what the
physical pixels are. Let's look at only a single dimension (x): say that the
minimum offset is x_min, and that the maximum extent (defined as the offset plus
the size) is x_max. Then the physical pixels these correspond to range from
x_min to x_max, so the size of the output is x_max - x_min, and the offset of
the output is x_min.

Yes, I envision the pixel stacks being transferred into a vector which is used
with psVectorStats.

Your algorithm looks OK, apart from the weights (see above).

[Paul Price]

Closing subsequent to release of SDRS-08, ADD-07.

[gusciora@…]

Please elaborate more on the Poisson weights. How are the weights to be
calculated (what is the formula or algorithm)? What should be done with them?

From the SDRS: "If zero is non-NULL and applyZeroScale is false, then the values
shall only be used in calculating the noise.". What does that mean? How should
the values be used in calculating the noise (hat is the formula or algorithm)?
What should be done with the noise once it is calculated.

A similar question applies to the scale psVector.

[Paul Price]

The formula for Poisson weights is provided above (comment 3), and I have
included it in the SDRS.

For "noise", read "Poisson variance".

I have updated the SDRS:

If the \code{zero} vector is non-\code{NULL} and \code{applyZeroScale}
is \code{true}, then the appropriate values shall be added to the
\code{inputs} before rejection is performed. If \code{zero} is
non-\code{NULL} and \code{applyZeroScale} is false, then the values
shall only be used in calculating the Poisson variances.

....

If the \code{gain} and \code{readnoise} are positive and non-negative
(respectively), then these shall be used to provide weights for the
combination using Poisson statistics: the Poisson variance is
$\sigma^{2 = (gf + r}2)/g$, where $\sigma^{2$ is the variance, $g$ is
the gain, $f$ is the value of the pixel, and $r$ is the read noise.}

[gusciora@…]

Let me try a different approach.

The module takes the pixels from N images and builds a "tmp" psVector. It adds
the elements from the zero psVector to the associated element in the "tmp"
vector. It multiplies the element in this "tmp" vector by the associated
element in the scale vector. (Incidentally, the SDR does not specify whether
the add or the multiply comes first. Even though they produce very different
results, both choices can be deemed "correct" from the requirements point of
view). It then calls psVectorStats on that tmp vector and sets the output pixel
to the stat value.

I assume the Poisson weight stuff means that this weight multiplies the input
pixel value. If so, when? Before the zero and scale vectors are applied, or
after? What is the formula for the weight? I assume it depends on that sigma^{2
formula you provided, but the description is not complete. Are you saying that
the Poisson variance *is* the weight?}

[Paul Price]

I've spelled it out more explicitly (see below). Please let me know if there is
still confusion.

[....]

If the \code{gain} and \code{readnoise} are positive and non-negative
(respectively), then these shall be used to provide weights for the
combination using Poisson statistics ($\sigma_i$ below).

In summary, pixels corresponding to the same physical pixel are
combined, having values $x_i \pm \sigma_i$. In the case that
\code{applyZeroScale} is \code{true}, then:
\begin{eqnarray}
x_i & = & s_i f_i + z_i
\sigma_i & = & (g x_i + r²⁾{1/2} / g
\end{eqnarray}
Where $f_i$ is the value of the pixel in image $i$, $s_i$ is the scale
applied to image $i$, $z_i$ is the zero offset applied to image $i$,
$g$ is the gain, and $r$ is the read noise. If scales are not
provided, they are set to unity; if zero offsets are not provided,
they are set to zero.

If \code{applyZeroScale} is \code{false}, then the values are:
\begin{eqnarray}
x_i & = & f_i
\sigma_i & = (g (s_i f_i + z_i) + r²⁾{1/2} / g
\end{eqnarray}
where the same symbols are used as above.

[gusciora@…]

How do I combine pixels having values $x_i \pm \sigma_i$?

Your explanation pertains to the case where "gain and readnoise are positive and
non-negative (respectively)". How do I combine pixels when they are not?
Simple addition (with the zero and scale vectors applied)?

[gusciora@…]

I have an idea. Should I put the sigma_i values in an error vector and let
psVectorStats() handle the errors when combining the pixels? Is that what you
had in mind?

[Paul Price]

Exactly. (:

[gusciora@…]

Yes! Okay, you can close this call.