I have a function like this $f:\mathbb{R}^2\to\mathbb{R}_+$
$$f(x, y) = e^{-\frac{x^2+y^2}{2}}$$
Its level sets $f(x, y) = c$ are simply circles centered at the origin $$ x^2 + y^2 = \log\left(\frac{1}{c^2}\right) \qquad \forall\, c\in\mathbb{R}_+ $$
Here's my question:
Usually we would take the integral of this function as a double integral where one dimension of integration is $x$ and one is $y$ $$ \int_{\mathbb{R}^2} f(x, y) dxdy = \int_{\mathbb{R}} \int_{\mathbb{R}} f(x, y) dx dy $$ Can we rewrite this integral as an integral over its level sets? For instance something like this $$ \int_{\mathbb{R}^2} f(x, y) dxdy = \int_{\mathbb{R}_+} \int_{x^2 + y^2 = \log\left(\frac{1}{c^2}\right)} f(x, y) dt dc $$
I am not sure how to write it and I am pretty sure the "inner" integral should be modified so that we are integrating in a "curvilinear" way along the circle.

In a way, it's something called a Lebesgue integral. Assuming that $f(x,y) > 0$ we have
$$ \int_{\mathbb{R}^2} f(x,y) dx dy = \int_0^\infty \mu\big(\{(x,y): f(x,y) > t\}\big) dt $$
where $\mu(A)$ is generally called the Lebesgue measure of a set, and in the case of 'nice' sets it's just the area; in this case, the area of the region where $f(x,y) > t$.
With the function $f(x,y) = e^{-\frac{x^2+y^2}{2}}$ we have: \begin{align}\{(x,y): f(x,y) > t\} &= \{(x,y): e^{-\frac{x^2+y^2}{2}} > t\} = \{(x,y): x^2+y^2 < - 2 \ln t\} =\\ &= \left\{\begin{array}{ll} \varnothing & \text{for }t\ge 1 \\ \{(x,y): x^2+y^2 < 2 |\ln t|\} & \text{for }0<t< 1\end{array} \right. \end{align} so $$ \mu\big(\{(x,y): f(x,y) > t\}\big) = \left\{\begin{array}{ll} 0 & \text{for }t\ge 1 \\ 2\pi| \ln t| & \text{for }0<t< 1\end{array} \right.$$ and you get
$$ \int_{\mathbb{R}^2} f(x,y) dx dy = \int_0^1 2\pi| \ln t| dt = 2\pi(t-t\ln t)\big|_0^1 =2\pi$$
The Lebesgue integral can be defined on other sets, not only $\mathbb R^n$, as long as you can measure the size of sets somehow.
ANOTHER APPROACH:
Another way you can try to calculate this integral is to parametrize $\mathbb R^2$ with new coordinates $x= x(c,t)$, $y=y(c,t)$ such that $$ f(x(c,t),y(c,t)) = c$$ Then we have $$ \int_{\mathbb{R}^2} f(x,y) dx dy = \iint f(x(c,t),y(c,t)) |J(c,t)| dc dt = \iint c |J(c,t)| dc dt$$ where $J(c,t)$ is the jacobian: $$ J(c,t) = \det\left[\begin{array}{cc} \frac{\partial x}{\partial c}& \frac{\partial x}{\partial t}\\ \frac{\partial y}{\partial c}& \frac{\partial y}{\partial t}\end{array}\right]$$
In this case, we could use coordinates $(c,t)$ such that $c\in(0,1]$, $t\in[0,2\pi]$ and $$ x= \sqrt{-2\ln c} \cos t$$ $$ y= \sqrt{-2\ln c} \sin t$$ Then $$ |J(c,t)| = \frac{1}{c} $$ so $$ \int_{\mathbb{R}^2} e^{-\frac{x^2+y^2}{2}} dx dy = \iint_{(c,t)\in(0,1]\times[0,2\pi]} c \frac{1}{c} dc dt = \int_0^{2\pi} \Big(\int_0^1 dc\Big) dt = 2\pi $$ In this case, the integral turns out to be very simple, but in general it's difficult to find parametrization $x= x(c,t)$, $y=y(c,t)$. This case was special because the level sets were concentric circles, so just a slight modification of well-known radial coordinates did the job.