We present an unconstrained tree tensor network approach to the study of lattice gauge theories in two spatial dimensions showing how to perform numerical simulations of theories in presence of fermionic matter and four-body magnetic terms, at zero and finite density, with periodic and open boundary conditions. We exploit the quantum link representation of the gauge fields and demonstrate that a fermionic rishon representation of the quantum links allows us to efficiently handle the fermionic matter while finite densities are naturally enclosed in the tensor network description. We explicit perform calculations for quantum electrodynamics in the spin-one quantum link representation on lattice sizes of up to 16x16 sites, detecting and characterizing different quantum regimes. In particular, at finite density, we detect signatures of a phase separation as a function of the bare mass values at different filling densities. The presented approach can be extended straightforwardly to three spatial dimensions.