Thanks. I've re-written the function to minimize the amount of copying 
(i.e. slicing) that is required. But now, I'm befuddled as to how to 
parallelize this function using Julia. Any suggestions?

Alex

function expensive_hat(S::Array{Complex{Float64},2}, mx::Array{Int64,2}, 
my::Array{Int64,2})

    samples = 64
        A = zeros(size(mx));
    B = zeros(size(mx));
    C = zeros(size(mx));

    for i = 1:samples
        Si = S[:,i];
        Sx = Si[mx];
        Sy = Si[my];
        Sxy = Si[mx+my];
        Sxyc = conj(Sxy);

                A +=  abs2(Sy .* Sx);
        B += abs2(sqrt(Sxyc .* Sxy));
        C += Sxyc .* Sy .* Sx;
    end

        ans = (A .* B ./ samples ./ samples, C./samples)
    return ans

end

data = rand(24000,64);
limit = 2000;

ix = int64([1:limit/2]);
iy = ix[1:end/2];
mg = zeros(Int64,length(iy),length(ix));
mx = broadcast(+,ix',mg);
my = broadcast(+,iy,mg);
S = rfft(data,1)./24000;

@time (AB, C) = expensive_hat(S,mx,my);
 

Reply via email to