function [blurIm] = blurImageVectors(testIm, testImSz, sigma)
  % [blurIm] = blurImageVectors(testIm, testImSz, sigma)
  % Input:
  %  testIm - N x K matrix with each of K images stored in the columns.
  %           N = nx*ny where the images have ny rows and nx columns. 
  %  testImSz = [ny, nx]
  %  sigma - stdev of gaussian blur. Here 6*sigma + 1 < min(nx,ny)
  
  sigmaSqr = sigma*sigma;
  gFiltSize = 2 * round(3.0 * sigma) + 1;
  x = [1:gFiltSize] - round((gFiltSize+1)/2);
  gFilt = exp(- x .* x / (2.0*sigmaSqr));
  gFilt = gFilt/ sum(gFilt(:));

  testIm = reshape(testIm, [testImSz, size(testIm,2)]);
  for k = 1:size(testIm,3)
    im = testIm(:, :, k);
    bIm = rconv2sep(im, gFilt, gFilt);
    blurIm(:,:,k) = bIm;
  end
  blurIm = reshape(blurIm, [prod(testImSz), size(testIm,3)]);
  
  return;