Hi Roland- I'm glad that helped. You're right about the iterative mean/std since the mask varies based on the current threshold which means you need to recalculate the entire thing anyway.
Using statsover does calculate all and then discard the ones you don't use. As far as a decision to optimize further, I suggest using the Benchmark module or Devel::NYTProf to measure the performance directly. Regards, Chris On Wed, Mar 12, 2014 at 10:33 AM, Schregle Roland HSLU T&A <[email protected]> wrote: > On Tue, 11 Mar 2014 23:28:20 +0100, Chris Marshall <[email protected]> > wrote: > >> Hi Roland- >> >> It looks like you've already gotten everything working but >> the final selective average. I've put comments inline with >> your text below on ways I would do things for efficiency if >> performance (speed or memory) were important. > > > Hi Chris, > > megathanks up front for your feedback -- you just catapulted our project > forward not insignificantly. ;^) > > >> Since glue() and append() all reallocate the entire pdl object to >> add one entry, it is more efficient to allocate by bigger chunks >> to reduce the overhead. Then just use .= slice assignment to >> add new frame data to the stack. > > > The images are generated in p parallel threads, so I could allocate a > [3,M,N,p]-dim PDL and drop them in there as the threads complete, then > glue/append them to the stack in one go before averaging. > > >> An array slice assignment clarifies the LHS to capture the desired >> 2 outputs from statsover: >> >> ($mean,$sigma) = ( $imgStackLum->mv(-1,0)->statsover )[0,6]; > > > How efficient is it to use statsover() here, considering I just need the > mean and stddev and ignore the remaining output? Does it actually skip > calculating the unwanted stats or are they simply thrown away? Would it be > faster to get mean and stddev manually? > > >> You could use bad values to perform the computation but another >> approach is to calculated the average by hand. The sum of the values >> is: >> >> $imgTotal = ( $imgStack * $mask(*3) )->mv(-1,0)->sumover; >> >> and the number of values is just the same sum but just over the >> mask with one subtlety: where the N==1 you need to avoid a >> division: >> >> $imgNum = $mask->mv(-1,0)->sumover; >> >> And the resulting average over the participating pixels is then: >> >> $imgComb = $imgTotal / $imgNum->lclip(1); > > > This is great -- and it works! I was tempted to use average() after > multiplying with $mask, until I realised it uses uniform rather than > per-pixel weights; needless to say, the results were even noisier than > without sigma clip! ;^) > > >> The usual things for iterative calculation of statistics is to keep >> the raw sum and sum-of-squares totals so you can update the >> values for each pixel in O(1) time rather than O(n+1). > > > Updating a raw sum for the mean is straightforward, but I can't see a way > around (mostly) redoing the stddev every time. Similarly, $imgTotal and > $imgNum depend on $mask, which also changes with every iteration (albeit > less frequently as more images are added). > > >> I'm assuming (hoping) that you are using a current version of PDL >> but the above should be correct for the last several releases. >> Maybe someone else will post with the example using bad values. > > > I'm using PDL v2.007 with Perl 5.10 and updated all the relevant modules via > CPAN. Your suggestions work fine! > > Craig posted an example using bad values; will look into that next. > > Many thanks again, much obliged! > > 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
