Speeding up this image smoothing algorithm

I have an image smoothing algorithm that is killing performance. It seems that there should be a better way to do this.

    win = 1-d profile.
    f = image of doubles.
    In = -2:2;
    lm = -2:2;  
    for i=3,n-2
        for j=3,m-2
           g(i,j) = sum(sort(reshape(f(i+In,j+lm),1,25)).*win);
        end
    end

... roughly speaking for each pixel we extract the 5X5 neighborhood, put it into a 1-d array, sort it and form the dot product with another vector.

At the image size I have, a median filter takes 7 seconds. This thing takes about 30. Incidentally, if I put it all into a gpu, it takes far longer. If I want to cling to this smoothing, how can I improve performance? Is there a clever way to vectorize this.

ideas welcome