Improving three-dimensional (3D) localization precision is of paramount importance for super-resolution imaging. By properly engineering the point spread function (PSF), such as utilizing Laguerre–Gaussian (LG) modes and their superposition, the ultimate limits of 3D localization precision can be enhanced. However, achieving these limits is challenging, as it often involves complicated detection strategies and practical limitations. In this work, we rigorously derive the ultimate 3D localization limits of LG modes and their superposition, specifically rotation modes, in the multi-parameter estimation framework. Our findings reveal that a significant portion of the information required for achieving 3D super-localization of LG modes can be obtained through feasible intensity detection. Moreover, the 3D ultimate precision can be achieved when the azimuthal index l is zero. To provide a proof-of-principle demonstration, we develop an iterative maximum likelihood estimation (MLE) algorithm that converges to the 3D position of a point source, considering the pixelation and detector noise. The experimental implementation exhibits an improvement of up to two-fold in lateral localization precision and up to twenty-fold in axial localization precision when using LG modes compared to Gaussian mode. We also showcase the superior axial localization capability of the rotation mode within the near-focus region, effectively overcoming the limitations encountered by single LG modes. Notably, in the presence of realistic aberration, the algorithm robustly achieves the Cramér-Rao lower bound. Our findings provide valuable insights for evaluating and optimizing the achievable 3D localization precision, which will facilitate the advancements in super-resolution microscopy.