Objective: Large-animal single-photon emission computed tomography (SPECT) systems are crucial for preclinical cardiovascular research. An accurate system matrix (or system response matrix) is essential for high-quality iterative image reconstruction. However, directly measuring the system matrix on large-animal SPECT systems is often prohibitively time-consuming due to the extensive field of view and high spatial resolution needs. To address the lengthy measurement process for the system matrix in large-animal SPECT, this work proposes and validates a calculation method based on two-dimensional (2D) Gaussian fitting. This method leverages the inherent continuity of the projection probability density function (PPDF) in the image domain. Instead of measuring the system response for each voxel, the proposed approach demonstrates that an accurate system matrix can be built by acquiring only sparse point-source measurement data. Using 2D Gaussian fitting, the method effectively models the system's spatial response and connects sparse data points to synthesize the full matrix. Methods: The complete dataset of point-source projections for 37 210 voxels was collected on the large-animal SPECT system, with a total acquisition time of 124 h. Down-sampled subsets at ratios of 1/4 and 1/9 were created from the full dataset to simulate accelerated protocols and to thoroughly evaluate the feasibility and robustness of the proposed method. The accuracy of the fitted PPDFs was quantitatively assessed against the fully measured ground truth using two metrics: the relative root mean square error (RRMSE) and the structural similarity index measure (SSIM). Additionally, to evaluate the reconstruction performance of the system matrix derived from the fitted PPDFs both qualitatively and quantitatively, 3.5 and 4 mm hot-rod phantom images were reconstructed. The performance of the proposed method was compared with several traditional approaches, including fully sampled direct measurement, barycentric Lagrange interpolation, and cubic spline interpolation. Results: Quantitative analyses showed exceptional fidelity in the estimated system matrices. The matrix computed under the 1/4 sparse sampling condition achieved an RRMSE of 0.023 ± 0.068 and an SSIM of 0.997 ± 0.003. Even with the more aggressive 1/9 sparse sampling, the method produced an RRMSE of 0.035 ± 0.095 and an SSIM of 0.996 ± 0.006. The system matrix generated with this method successfully resolved 3.5-mm hot rods. Under the 1/4 sparse sampling, the reconstructed images of 3.5 and 4 mm rods had RRMSE values of 0.795 and 0.654 against the ground truth, and SSIM values of 0.981 and 0.988. For the same sampling, the reconstructed images showed RRMSE values of 1.042 and 0.797, with SSIM values of 0.971 and 0.981. These imaging results were visually and quantitatively superior to those obtained with other methods. Importantly, the computational time remained efficient, with calculations taking only 1.9 hours and 1.8 hours for the 1/4 and 1/9 sparse datasets, respectively. Conclusions: The proposed 2D Gaussian fitting approach effectively overcomes the traditional limitations in generating system matrices for large-animal SPECT systems. It significantly reduces measurement time and computational costs without sacrificing tomographic image quality. This method presents a practical and efficient solution for acquiring precise system matrices in large-animal SPECT imaging.