CODEWITHSUNDEEP
julia
aramirezreyes
aramirezreyes

Reputation: 1365

Julia equivalent to movmean(array,window,dims)?

I am using julia 0.7.0, NCDatasets.jl and Images.jl on a linux box to analize a dataset of around 80GB. I don't load a lot of variables and the first step is to do the equivalent of matlab's

a = moveman(movemean(movemean(array,window,1),window,2),window,4))

where array is a (256,256,80,600) float array. For this, I am trying the line: 

filtered = imfilter(array, centered(ones(window_h,window_h,1,window_t)/(window_t*window_h*window_h)),Inner())

However, this results in terabytes of allocations, which ends up using all my memory and taking ages. The matlab line works just fine and use an insignificant amount of time compared to that of my julia line, which suggests I am doing something in a non-optimal way.

Could someone provide any insight?

Upvotes: 2

Views: 608

Answers (2)

aramirezreyes
aramirezreyes

Reputation: 1365

To answer my own question, based on the discussion at Julia discourse: I continued using the Images package, in particular ImageFiltering in the following way. First I define the kernel used to smooth. This kernel will be used by computing the correlation between it and the array that we are filtering.

The difference by using factores kernels is that the each filter will be applied separately, which changes the number of operations from

window_h x window_h x window_t

to

window_h + window_h + window_t

As explained in the docs.

Note that the kernel uses [1.0] in the third dimension, because my array is a 4-dim array and I am smoothing on the first two and the fourth dimension.

using ImageFiltering  

function kernel4d_2(window_h,window_t)
    kernel_h = ones(window_h)/window_h
    kernel_t = ones(window_t)/window_t
    return kernelfactors((kernel_h, kernel_h, [1.0], kernel_t))
end

Then I defined a function to apply this kernel as a filter and return the filtered array.

function filter_array(array,window_x,window_t)
        filtered = imfilter(array, kernel4d_2(window_x,window_t))
end

This allows to filter the array as:

filtered = filter_array(unfiltered,window_x,window_t)

Upvotes: 1

Jason D
Jason D

Reputation: 303

not quite familiar with matlab, guess it's moving average?

then it's linear, and to do movemean(movemean(movemean...

you may calculate an equation instead, like

( 3*array[current] + 3*array[current-1] + 2*array[current-2] )/8

and go through the array

Upvotes: 1

Related Questions