Hyperspectral unmixing is concerned with how to learn the constituent materials (called endmembers) and their corresponding proportions (called abundances) from mixed pixels. Nonnegative tensor factorization (NTF) is powerful in addressing the unmixing problem owing to the ability to preserve spatial inherent information. However, the neglect of local spatial information, nonlinear mixtures, and spectral variability degrades the quality of unmixing result. In this article, we propose a superpixel-based low-rank tensor factorization (SLRTF) for blind nonlinear hyperspectral unmixing. Our method first generates numerous superpixels and processes them into regular cubes to exploit the local spatial information. Then, the resulting superpixel cubes are jointly decomposed using a nonlinear low-rank tensor factorization based on a generalized bilinear model (GBM) to account for the second-order scattering phenomenon. Notably, each superpixel cube yields an endmember matrix to consider spectral variability. Furthermore, to make full use of the consistent and complementary spectral information, a double-weighted endmember (DWE) constraint is designed for adaptively learning the consensus endmember matrix. The tensor nuclear norm constraint is enforced on the abundance tensor to characterize high spatial correlation. Finally, the proposed method is compared with several state-of-the-art methods on synthetic and real datasets. Experimental results demonstrate that the proposed method offers superior performance for hyperspectral unmixing.