This study employs the two-dimensional Legendre polynomial fitting (2-D LPF) method to fit M2 tidal harmonic constants from satellite altimetry data within the region of 53°E–131°E, 34°S–6°N, extracting internal tide signals acting on the sea surface. The M2 tidal harmonic constants are derived from the sea surface height (SSH) data of the TOPEX/Poseidon (T/P), Jason-1, Jason-2, and Jason-3 satellites via t-tide analysis. We validate the 2-D LPF method against the 300 km moving average (300 km smooth) method and the one-dimensional Legendre polynomial fitting (1-D LPF) method. Through cross-validation across 42 orbits, the optimal polynomial orders are determined to be seven for 1-D LPF, and eight and seven for the longitudinal and latitudinal directions in 2-D LPF, respectively. The 2-D LPF method demonstrated superior spatial continuity and smoothness of internal tide signals. Further single-orbit correlation analysis confirmed generally higher correlation with topographic and density perturbations (correlation coefficients: 0.502, 0.620, 0.245; 0.420, 0.273, −0.101), underscoring its accuracy. Overall, the 2-D LPF method can use all regional data points, overcoming the limitations of single-orbit approaches and proving its effectiveness in extracting internal tide signals acting on the sea surface.