## 1 Introduction  
Marine seismic data are accompanied by strong multiples, and addressing these multiples has become a crucial aspect of seismic processing. Traditional migration algorithms focus on primaries, and multiples are treated as noise and are removed from the data before seismic migration. However, in contrast to other interference noise, multiples represent true seismic responses that have experienced multiple scattering and contain information about subsurface structures. Multiples can be integrated into seismic migration and inversion as valuable signals. Compared to primaries, surface-related multiples can be regarded as plane sources at the free surface, offering an extended lateral illumination range relative to the sparser point sources. In terms of reflection angles and propagation paths, surface-related multiples have smaller reflection angles and longer propagation paths. Higher-order multiples generate even smaller reflection angles, leading to increased vertical resolution of the imaging results. Internal multiples are produced at strongly reflective interfaces in the subsurface, such as salt domes, and can be interpreted as additional secondary scattering sources between layers. This provides supplementary illumination beneath substantial shielding layers, such as salt domes. Many studies have demonstrated the advantages of utilizing multiple imaging [1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11, 12, 13, 14]. Seismic migration imaging technology has played a positive role in obtaining underground reflection coefficients and geological structure investigations [15; 16]. However, the migration results represent only approximate values of the subsurface reflection coefficients. To further seek the accurate distribution of underground structures, the least-squares technique is widely applied in migration imaging and tomographic inversion [17; 18; 19; 20; 21], yielding excellent results. LSM typically considers scenarios involving only primary reflections. However, when data contain both primary and multiple reflections, multiples generate a substantial amount of crosstalk artifacts, which pose significant constraints on the practical application of imaging techniques. Coherent artifacts of surface-related multiples primarily arise from improper correlation imaging between unrelated multi-order wavefields. Internal multiples have complex and variable paths, complicating the simulation of internal multiple paths during forward propagation and posing challenges in achieving correct imaging. Currently, various methods are utilized to suppress the crosstalk. One category of methods focuses on migration imaging conditions, such as migration using deconvolution-type imaging conditions [22; 23; 24] and stereographic imaging conditions [25]. Deconvolution is effective for short-period crosstalk, but it underperforms in suppressing strong coherent crosstalk events [26]. Another category of methods includes least-squares migration methods, which can be further classified into those that separate multiples for least-squares migration and those that utilize the simultaneous imaging of joint primaries and multiples for least-squares migration. Specifically, in least-squares reverse time migration of multiples (L8 RTMM) containing primaries and surface-related multiples [27; 28], surface-related multiples serve as plane sources, which, together with point sources, form the total downward wavefield for migration. Subsequently, the least-squares strategy is employed to suppress the crosstalk. The generalized internal multiple imaging (GIII) method can simulate high-order internal scattering through interference and can subsequently develop least-squares suppression of artifacts [29; 30]. Full-wavefield migration (FWM) offers various strategies for combining primary and multiple imaging [31; 32; 33; 34; 35] and uses closed-loop one-way wave operators and inversion strategies to simultaneously process data containing surface-related multiples and internal multiples.