function color = pixelValue ( image , x, y, method )
s = size(image);
% pixel value at real coordinates
% do the interpolation
switch ( method )
    case 'nearest'
        % nearest neighbour interpolation
        fx = floor(x+0.5);
        fy = floor(y+0.5);
        if inImage(s, fx, fy)
            color = image(fy,fx);
        else
            color = 0;
        end
    case 'linear'
        % Do bilinear interpolation
        i = floor(y);
        j = floor(x);
        i1 = i+1; %ceil(y);
        j1 = j+1; %ceil(x);
        dx = x - j;
        dy = y - i;

        if inImage(s, j, i) && inImage(s, j1, i1)
            color = (1-dx) * (1-dy) * image(i, j);
            color = color + dy * (1 - dx)*image(i1,j);
            color = color + (1-dy) * dx * image(i,j1);
            color = color + dx * dy * image(i1,j1);
        else
            %return random noise, useful for debugging to distinguish
            %values out of bounds from 'real' black pixels of the original
            %color = rand;
            
            %but this looks more professional..:
            color = 0;
        end
end %end switch