You could generate a weighting variable and use the sum over the weighting
variable, instead of average(). But you're better off using BAD values, which
are awesome. You can set a given element to a special value that indicates
it's bad, and it will be ignored by all the statistical routines. Also, if
you're using glue() iteratively, you'll do marginally better to accumulate n of
your (3,N,M) images into a Perl array, and then PDL them all at once. Like this:
$imgStack = pdl( double, map { rpic($_) } @files ); # 3,N,M,n
$img_NMn = $imgStack->average; # N,M,n;
$img_NMn->badflag(1);
$img_nNM = $img_NMn->mv(2,0);
do {
my($mean, undef, undef, undef, undef, undef, $sigma) =
$img_nNM->statsover;
$dex = whichND( abs($img_NMn-$mean) > $sigma);
$img_NMn->indexND($dex) .= $img_NMn->badvalue;
} until($img_NMn->nelem==0);
On Mar 11, 2014, at 2:44 PM, Schregle Roland HSLU T&A <[email protected]>
wrote:
> Hi all,
>
> I'm new to PDL and semi-proficient with Perl at best. I'm posting this a last
> resort after *lots* of frustration while trying to combine a set of RGB
> images of dimensionality [3, N, M] using a selective averaging technique
> often used in astrophotography called sigma clip (see
> http://www.ccdware.com/support/presentations/advanced_image_combine_gralak20081116.pdf).
>
> The images are generated iteratively and the set size is not known in advance.
>
> I append the images into a 4D array on every iteration using glue(). Of
> course I'm open to alternatives...
>
> I'm seriously stuck on actually combining the [3, N, M, n] image stack after
> n iterations using sigma clipping. This requires getting the [N, M] mean and
> standard deviation (sigma) from the pixel luminance (average over RGB), then
> averaging the stacked images using only those pixels whose deviation from the
> mean is within some factor of sigma. This effectively rejects outliers.
>
> Here's what I have so far...
>
> ----------------------------------------------------------------------------------------------------
> use constant SIGMACLIP => 0.5;
> my $imgStack;
>
> do {
> ...
> # Append new image $img to stack, adding 4th dimension
> $imgStack = defined $imgStack ? $imgStack -> glue(3, $img) : $img;
> ...
>
> # Get mean and std deviation for RGB average (= luminance) of all images
> my $imgStackLum = average($imgStack);
> # [N, M, n] --> [n, N, M]
> my ($mean, undef, undef, undef, undef, undef, $sigma) =
> statsover($imgStackLum -> reorder(2,0,1));
>
> # Add 3rd dim to mean for subtraction
> my $imgStackDiff = $imgStackLum - $mean -> dummy(2);
>
> # Get mask array containing non-clipped pixels
> $sigma *= SIGMACLIP;
> my $mask = (abs($imgStackDiff) < $sigma);
>
> # Selectively combine $imgStack into selectively averaged [3, N, M] image
> based on $mask
> $hdrComb = where($hdrStack, $mask)... -> average; ???
>
> ...
> } until (converged);
> ----------------------------------------------------------------------------------------------------
>
> Note that I have to redo this for every iteration, which is anything but
> efficient. Again, I'm open to suggestions for optimisation... once it's
> working.
>
> This is all pretty standard stuff, so nothing exotic, and I'll bet somebody's
> done it in PDL. I couldn't find any hints in the archives tho. (Apparently
> the search function in the recent archives is broken, btw).
>
> Any help appreciated. Many thanks & best regards,
>
> --Roland
>
>
> --
> Dr. Roland Schregle
> Senior Research Associate
>
> T direct: +41 41 349 39 77
> [email protected]
>
> Lucerne University of Applied Sciences and Arts
> School of Engineering and Architecture
> CC Envelopes and Solar Energy (EASE)
> Technikumstrasse 21, CH-6048 Horw
> T +41 41 349 33 11, F +41 41 349 39 60
> www.hslu.ch/ccease
>
> _______________________________________________
> Perldl mailing list
> [email protected]
> http://mailman.jach.hawaii.edu/mailman/listinfo/perldl
>
_______________________________________________
Perldl mailing list
[email protected]
http://mailman.jach.hawaii.edu/mailman/listinfo/perldl
