Chris nailed this, as always. You're probably running out of memory with the range() call, simply because range() is implemented in the most boneheaded way possible, which is fine for smaller data sets -- but for a 1k x 1k image you're already looking at a gigapixel output from range. Chris's initiative to remove the 2GB Perl memory restriction should fix that, but meanwhile you're in territory where memory use actually matters -- and therefore out of the domain in which range() is meant to operate.
I suggest using explicit looping inside an inline PP routine -- that way you can keep everything in cache and get optimal performance. Here is an example inline PP routine that also uses C helper function -- it helps to have a template to start from, but it's pretty straightforward to do. (The only hassle with Inline::PP is getting the darned thing to compile -- you have to drill down into the temp directory to find any error messages made by the compiler).
hilbert.pdl
Description: Binary data
On Sep 25, 2011, at 7:36 AM, chm wrote:
> On 9/25/2011 8:40 AM, Jarrod Miller wrote:
>> Hi folks.
>> I have a very simple question about threading that I'm sure won't be
>> difficult for a seasoned piddle twiddler..
>>
>> I have a 2D matrix which is essentially a greyscale image.
>> I need to implement a moving window over the image to perform a few
>> transformations based on surrounding pixel values. For now I'll likely just
>> be averaging the pixel intensity around each pixel so I can build a
>> threshold.
>> At the moment I'm using this code to get my moving window(s), inspired by
>> the threads page on the wiki: $p->range($p->ndcoords()-15, 31,
>> 'm')->reorder(2,3,0,1); This works fantastic -- except for the part where I
>> run out of RAM and perl segfaults.
>
> range is very powerful but memory intensive. In your case,
> the piddle created is almost 1000X larger than the original.
> If perl runs out of memory then it exits. I'm not sure why
> you are seeing segfaults.
>
>> So I was wondering how I should go about doing this. My current "workaround"
>> is to use perl loops, and just move the window one slice at a time, adding
>> and subtracting the columns as I get to them, so I can calculate a
>> relatively fast average. It works, but I'm very sure this isn't the 'right'
>> way to do it..
>> Any help would be greatly appreciated.
>
> Looping over smaller regions/slices is one approach but if
> you are using range as above, you'll still see the 1000X
> memory increase for the operations.
>
> The simplest approach is to use one of the existing
> PDL routines for 2D data. I.e.,
>
>
>> $ pdldoc -a 2d
>> PDL::Graphics2D Module: An object oriented interface to PDL graphics
>> PDL::Image2D Module: Miscellaneous 2D image processing functions
>> applywarp2d Transform a set of points using a 2-D polynomial mapping
>> bilin2d Bilinearly maps the first piddle in the second. The
>> interpolated values are actually added to the second piddle which is
>> supposed to be larger than the first one.
>> box2d fast 2D boxcar average
>> cc8compt Connected 8-component labeling of a binary image.
>> centroid2d Refine a list of object positions in 2D image by centroiding
>> in a box
>> conv2d 2D convolution of an array with a kernel (smoothing)
>> fitwarp2d Find the best-fit 2D polynomial to describe a coordinate
>> transformation.
>> hi2d Plot image as 2d histogram (not very good IMHO...)
>> histogram2d Calculates a 2d histogram.
>> imag2d Display a 2-D image in a figure window
>> imagrgb 2D RGB image plot (see also imag2d)
>> imagrgb3d 2D RGB image plot as an object inside a 3D space
>> index `index' and `index2d' provide rudimentary index indirection.
>> index2d `index' and `index2d' provide rudimentary index indirection.
>> inner2d Inner product over 2 dimensions.
>> max2d_ind Return value/position of maximum value in 2D image
>> med2d 2D median-convolution of an array with a kernel (smoothing)
>> med2df 2D median-convolution of an array in a pxq window (smoothing)
>> patch2d patch bad pixels out of 2D images using a mask
>> patchbad2d patch bad pixels out of 2D images containing bad values
>> polyfill fill the area inside the given polygon with a given colour
>> polyfillv return the (dataflown) area of an image within a polygon
>> polyfit Fit discrete data in a least squares sense by polynomials in
>> one variable. Handles threading correctly--one can pass in a
>> 2D PDL (as `$y') and it will pass back a 2D PDL, the rows of
>> which are the polynomial regression results (in `$r'
>> corresponding to the rows of $y.
>> rescale2d The first piddle is rescaled to the dimensions of the second
>> (expanding or meaning values as needed) and then added to it
>> in place. Nothing useful is returned.
>> rot2d rotate an image by given `angle'
>> twiddle Enable GUI interaction with a FreeGLUT display window.
>> warp2d Warp a 2D image given a polynomial describing the *reverse*
>> mapping.
>> warp2d_kernel Return the specified kernel, as used by warp2d
>> whistogram2d Calculates a 2d histogram from weighted data.
>
>
> conv2d seems appropriate for the current problem.
>
>> $ pdldoc conv2d
>> Module PDL::Image2D
>> conv2d
>> Signature: (a(m,n); kern(p,q); [o]b(m,n); int opt)
>>
>> 2D convolution of an array with a kernel (smoothing)
>>
>> For large kernels, using a FFT routine, such as fftconvolve() in
>> "PDL::FFT", will be quicker.
>>
>> $new = conv2d $old, $kernel, {OPTIONS}
>>
>> $smoothed = conv2d $image, ones(3,3), {Boundary => Reflect}
>>
>> Boundary - controls what values are assumed for the image when kernel
>> crosses its edge:
>> => Default - periodic boundary conditions
>> (i.e. wrap around axis)
>> => Reflect - reflect at boundary
>> => Truncate - truncate at boundary
>>
>> Unlike the FFT routines, conv2d is able to process bad values.
>
> You can also use PDL::PP code and Inline::Pdlpp to implement
> your own window computations directly for the threading engine.
>
> Cheers,
> Chris
>
> _______________________________________________
> 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
