When changing variables in multiple integrals one needs to find the Jacobian matrix of the transformation
$$\int dx_1\cdots dx_n f(x_1,...,x_n)=\int dy_1\cdots dy_n f(y_1,...,y_n)\left|\frac{\partial(x_1,...,x_n)}{\partial(y_1,...,y_n)}\right|$$
However, in the problem I am facing I don't understand why this does not work. I have to evaluate the integral
$$\int_0^1 dx_1\cdots dx_{2n} \delta(x_{i_1}-x_{i_1+1}+x_{i_2}-x_{i_2+1})\cdots \delta(x_{i_{n-1}}-x_{i_{n-1}+1}+x_{i_{n}}-x_{i_{n}+1})$$
with non-repetitive indices $i_1,...,i_n$. I want to substitute $y_i=x_i-x_{i+1}$ in the above integral. For this the Jacobin becomes simply $\frac{1}{2}$ and the multiple-integral factorizes to double integrals, but the result I am getting is definitely wrong. For simplicity I assume $i_{k-1}\pm1\neq i_k$ for any $k$. There must be some combinatorial problem of how the indices are arranged but I don't see it anywhere. Most likely I need something similar to Integration of Dirac delta over finite interval, but for arbitrary indices.
Any help for even small n<6 is helpful. Thanks in advance!