Calculating the number of solutions of integer linear equations Let $N$ be a natural number. Consider the following set of matrices whose entries are non-negative integers:
$$X_N:=\left\{(c_{ij})_{i,j=1}^4\in M_4(\mathbb{Z}_{\geq 0})\bigg| \sum_j c_{1j} = \sum_i c_{2j}, \sum_i c_{i1}=\sum_i c_{i2}, \sum_{i,j}c_{ij} = N\right\}.$$
Write $f(N) = |X_N|$. 
Question: Is there an algorithm / theorem which gives an explicit formula for $f(N)$ as a function of $N$? It should be polynomial. Without the constraints about the sum of the elements of the first two columns is equal and the sum of the elements of the first two rows is equal, this will just be $\binom{N+15}{15}$ which is a polynomial of degree 15. I am looking for a similar formula for the number of elements of $X_N$.
 A: OEIS entry A001496 has a formula and a generating function.
A: Let $r:=\sum_j c_{1j}$, $c:=\sum_i c_{i,1}$, and $s:=c_{11} + c_{12} + c_{21} + c_{22}$. Clearly, $s\leq N$, $2r\geq s$, $2c\geq s$, and $2r+2c-s \leq N$.
Next, let $u:=c_{11}$, $v:=c{12}$, $w:=c_{21}$, and so $c_{22}=s-u-v-w$. Then $u+v\leq r$, $w+(s-u-v-w)\leq r$, $u+w\leq c$, $v+(s-u-v-w)\leq c$, further implying $u+v\geq s-r$ and $u+w\geq s-c$.
Combining these inequalities together, we get
$$f(N) = \sum_{s=0}^N \sum_{r=\lceil s/2\rceil}^{ \lfloor N/2\rfloor} \sum_{c=\lceil s/2\rceil}^{\lfloor (N+s)/2\rfloor-r} \binom{N-(2r+2c-s)+3}{3} \sum_{u=0}^{\min\{ s,r,c\}} \sum_{v=\max\{0,s-r-u\}}^{\min\{ r-u, s-u, c\}} \sum_{w=\max\{0,s-c-u\}}^{\min\{ c-u, s-u-v, r\}} (1+r-u-v)(1+r+u+v-s)(1+c-u-w)(1+c+u+w-s),$$
where the binomial coefficients stands for the number of ways to fill up the lower-right $2\times 2$ minor (with the sum of elements equal $N-(2r+2c-s)$), and the product of terms in parentheses stands for number of suitable values for $c_{ij}$ for the upper-right and lower-left $2\times 2$ minors. 
