function exp = expandFunc(image, method, param1)

% generic expand function

if (strcmp(method,'burt'))
   W = param1;
   exp = upsample(image);
   exp = replicateConv(4*exp, W);
   exp = round(exp);
end

if (strcmp(method,'haar'))
   exp = upsample(image, 'duplicate');
end

if (strcmp(method,'ideal'))
   exp = 4*upsample(image);
   
   fft = fft2(exp);
   [row, col] = size(exp);
   
   xx = (-row/2+1):row/2;
   yy = (-col/2+1):col/2;
   [x, y] = meshgrid(xx,yy);
   mask = [and(abs(x)<row/4, abs(y)<row/4)];
%   mask = [(x.^2+y.^2) < row^2/16];
   mask = swapBlocks(mask);
   fft = fft.*mask;
   
   exp = real(ifft2(fft));
end

if (strcmp(method,'MMF'))
   exp = MIF(image);
end

if (strcmp(method, 'opt_lin'))
  exp = linear_interpolate(image);
end

if (strcmp(method, 'opt_cub'))
  exp = cubic_interpolate(image);
end

if (strcmp(method, 'spline'))
   order = param1;
   f1 = BsplineFilter(order);
   invf1 = invertFilter(f1, 512);
   u = binomialFilter(order);
   
   BSimage = replicateConv(image, invf1'*invf1);
   BSexp = upsample(BSimage);
   BSexp = replicateConv(BSexp, u'*u);
   exp = replicateConv(BSexp, f1'*f1);
end

if (strcmp(method, 'opt_NL'))
   coeffs = param1;
   part = partition(image);
   rimage1 = NLpredictor(image, coeffs(:,:,1), part);
   rimage2 = NLpredictor(image, coeffs(:,:,2), part);
   rimage3 = NLpredictor(image, coeffs(:,:,3), part);
   rimage4 = NLpredictor(image, coeffs(:,:,4), part);
   exp = zeros(2 * size(image));
   exp(1:2:size(exp,1), 1:2:size(exp,2)) = rimage1;
   exp(2:2:size(exp,1), 1:2:size(exp,2)) = rimage2;
   exp(1:2:size(exp,1), 2:2:size(exp,2)) = rimage3;
   exp(2:2:size(exp,1), 2:2:size(exp,2)) = rimage4;
end

if (strcmp(method, 'multi'))
   bestmethod = param1;
   nmethods = 10;
   [row col] = size(image);
   row = 2 * row;
   col = 2 * col;
   
   a = 0.4;
   w1 = [1/4-a/2; 1/4; a; 1/4; 1/4-a/2];
   W1 = w1*w1';
   a = 0.6;
   w2 = [1/4-a/2; 1/4; a; 1/4; 1/4-a/2];
   W2 = w2*w2';
   load NLcoeffs;
   
   rimage = zeros([row col nmethods]);
   rimage(:,:,1) = expandFunc(image, 'burt', W1); disp('burt  a=0.4');
   rimage(:,:,2) = expandFunc(image, 'burt', W2); disp('burt  a=0.6'); 
   rimage(:,:,3) = expandFunc(image, 'ideal'); disp('ideal');
   rimage(:,:,4) = expandFunc(image, 'haar'); disp('haar');
   rimage(:,:,5) = expandFunc(image, 'MMF', 3); disp('MMF');
   rimage(:,:,6) = expandFunc(image, 'spline', 1); disp('linear spline');
   rimage(:,:,7) = expandFunc(image, 'spline', 3); disp('cubic spline');
   rimage(:,:,8) = expandFunc(image, 'opt_lin', 15); disp('optimal linear');
   rimage(:,:,9) = expandFunc(image, 'opt_cub', 15); disp('optimal cubic');
   rimage(:,:,10) = expandFunc(image, 'opt_NL', NLcoeffs); disp('optimal non-linear');
   
   part = partition(image);
   meth = zeros([row col]);
   
   bm1 = bestmethod(1, :);
   bm2 = bestmethod(2, :);
   bm3 = bestmethod(3, :);
   bm4 = bestmethod(4, :);
   meth(1:2:row, 1:2:col) = bm1(part+1);
   meth(1:2:row, 2:2:col) = bm2(part+1);
   meth(2:2:row, 1:2:col) = bm3(part+1);
   meth(2:2:row, 2:2:col) = bm4(part+1);
   imshow(meth/10);
   
   part = reshape(part, [row*col/4 1]);
   
   exp = zeros([row*col 1]);
   meth = reshape(meth, [row*col 1]);
   
   for i = 1:nmethods
      rim = reshape(rimage(:,:,i), [row*col 1]);
      [tempo mask] = find([meth == i]');
      exp(mask) = rim(mask);
   end
   
   exp = reshape(exp, [row col]);
end