@yuggib, ok I have slightly modified the equation so that the convolution is well defined as it only act on common arguments x and y respectively. However I still encountered an issue.
Define 2D dirac delta, $\delta^2(x,y)=\delta(x)\delta(y)$. Then
$$(0e^{-y^2}+e^{-x^2})*\delta^2(x,y)=\int_{-\infty}^{\infty}\int_{-\infty}^{\infty}(0e^{-y'^2}+e^{-x'^2})\delta(x-x')\delta(y-y')dx'dy'$$
$$=\int_{-\infty}^{\infty}(0e^{-0^2}+e^{-x'^2})\delta(x-x')dx'=0e^{-0^2}+e^{-0^2}=1$$
$$\mathcal{F}(\delta^2(x,y))=1$$