K = A' * A;
[row_num,~] = size(f_noisy);
f = reshape(f, [row_num^2,1]);
f_noisy_tmp = reshape(f_noisy, [row_num^2,1]);
coeff = K \ (A' * f_noisy_tmp);
new_f = A * coeff;
scatter3(X(:),Y(:),new_f(:));

disp('error in relation to the original');
error = sum((f - new_f).^2);
disp(error);
disp('error in relation to the noise');
error = sum((f_noisy_tmp - new_f).^2);
disp(error);