Inertia tensor in non-cartesian coordinates given a rigid body $K$, I always had seen the formula
\begin{equation}
I_{ij} =  \int_K[\mathbf{x}^2\delta_{ij} - x_ix_j]\rho(\mathbf{x})\mathrm{d}^3\mathbf{x}
\end{equation}
for the inertia tensor. But in my derivation this seems to be correct only in cartesian coordinates (where $g_{ij} = \delta_{ij}$).
Here are my steps:
The angular momentum is
\begin{equation*}
        \mathbf{L} = \int_K (\mathbf{x} \times (\omega\times\mathbf{x})) \rho(\mathbf{x}) \mathrm{d}^3\mathbf{x}
\end{equation*}
written in components:
\begin{equation*}
        L^i = \int_K\left[\omega^i(x_jx^j) - x^i(x_j\omega^j)\right]\rho(\mathbf{x})\mathrm{d}^3\mathbf{x}
\end{equation*}
Now if I define:
\begin{equation*}
        I^i_{\ j} = \int_K[x_kx^k\delta^i_j - x^ix_j]\rho(\mathbf{x})\mathrm{d}^3\mathbf{x}
\end{equation*}
I get the relation:
\begin{equation}
L^i = I^i_{\ j}\omega^j
\end{equation}
Which is okay.
But if I write it by lowering the indices I get:
\begin{equation*}
        I_{ij} = g_{il}I^l_{\ j} = \int_K[x_kx^kg_{ij} - x_ix_j]\rho(\mathbf{x})\mathrm{d}^3\mathbf{x}
\end{equation*}
Which has a metric in place of the delta in the first formula seen in books.
I wonder if books give for granted that you work in cartesian coordinates .
Am I wrong?
EDIT:
I noticed that in my derivation I pushed $g_{ij}$ inside the integral while it actually may depend on coordinates. So I guess either one of $I^i_{\ j}$ or $I_{ij}$ I got is wrong? What is the correct definition?
 A: You are correct. The first formula does not work in general coordinates because $\delta_{ij}$ is not a tensor; it only transforms as one under orthogonal transformations. In general, it should be replaced by the metric.

Answering the rest of your question: First of all I missed something, which is that in general coordinates you should have $\sqrt{g}$ (where $g = \det g_{ij}$) inside the integral. But regardless of that, if we take the basic definition of the inertia tensor to be $L^i = I^i{}_j \omega^j$, we can then define $I_{ij} = g_{ik} I^k{}_j$, to be used in formulas like $E = 1/2 I_{ij} \omega^i \omega^j$. Therefore, the metric should go outside of the integral.
A: If you open any classical-mechanics text and it has $x^i,x_j,x,y,z$ or similar in itself, you should always assume that these are Cartesian coordinates and that the formulas won't work in any other coordinates system. 
Rewriting your formulas for a general coordinate system is easy when you are referring to components of vectors (pseudo-vectors) and in any expression you only refer to objects at a single location. Then the recipe is to replace every $\delta_{ij}$ with the curvi-linear metric $g_{ij}$, the Levi-Civita symbol as $\epsilon_{ijk} \to \epsilon_{ijk}\sqrt{\mathrm{det}(g)}$, and the vector components transformed by the components of the Jacobian of the transformation.
The case of the tensor of inertia is special, though. This is because it involves something which is called the position or distance "vector" $\vec{x}$. However, $\vec{x}$ does not transform as any other vector (with the Jacobian), it should be rather understood as a "convenient function triplet that kind of reminds us of a vector". This is because the distance "vector" only transforms as a vector only with respect to coordinate rotations. Consequently, $I^{ij}$ is not a tensor, it does not transform as one, and its indices should not be raised and lowered with the metric. If in doubt, ask yourself: At which point should one evaluate the metric components used to raise or lower the indices of $I^{ij}$? At the surface of the body you integrated over? At its center of mass? Obviously you should not do either. Once again, the moment of inertia only transforms as a tensor with respect to rotations of the Cartesian coordinate system.
This is not to say that the moment of inertia is always completely horrible in the context of curvilinear coordinates. This is because its real physical meaning is that it refers to inertia about a rotation axis, and once you have fixed which axis interests you (independent of the coordinate system), you can transform the corresponding integral as a usual scalar integral. For instance, the $I^{zz}$ component gives you the moment of inertia for rotation about the $z$-axis. In polar coordinates $R = \sqrt{x^2+y^2},\phi = \arctan(x/y),z=z$ the original $I^{zz}$ component can then be computed as
$$I^{zz} = \int_{0}^{2\pi}\!\!\!\!\int_0^\infty\!\!\!\int_{-\infty}^\infty\!\!\!\! R^3 \rho(R,\phi,z)dR d\phi d z  $$
This integral can be much simpler if, for instance, the body is rotationally symmetric about $z$ (the density is independent of $\phi$) etc.
