#Reviewer1

This is a review for the manuscript entitled “The crustal magma reservoir geometry and seismic activity beneath the Wudalianchi volcano in northeast China: Implications for the multilevel magmatic system” by Lü et al. This study conducts a new high-resolution model using a finite-frequency, full-wave method with a dense seismic network. The results show clear improvements in imaging the crustal magmatic structure of the Wudalianchi volcanic region compared to previous images. Notably, they identify a distinct magma reservoir in the middle-lower crust and propose potential connections within the magmatic system, from its deeper root to the shallow pre-eruptive reservoir. The data, methods and results appear robust, and the manuscript is interesting and well written. I have some suggestions and therefore suggest a minor revision. Please see the detailed comments below:

In the method section, Vp in the inversion helps with constraints of both shallow and deep structure. I wonder what the shallow Vp structure looks like. Can the author comment on this?

R. Thanks for your suggestion. Rayleigh waves are most sensitive to P-wave velocity at shallow depths, and a strong correlation is observed between the P-wave velocity anomaly and the regional faults at 3 km depth. However, the P-wave velocity structure is not well constrained at 11 km depth in comparison with the shear-wave structure due to the limited resolution. Text and figures have been updated in the revised manuscript.

Currently, the melt fraction of the middle-lower crustal magma chamber is calculated based on velocity perturbation instead of absolute velocity. Based on a dense seismic network, this study conducts a higher resolution model, and the full-wave method should provide a much better constraints to the absolute velocity. The result in absolute velocity may give a more precise estimate of the melt fraction

R. Thanks for your suggestion. We have added the melt fraction estimation based on absolute velocity. Figure has been updated in the revised manuscript.

In Figure S1, the temporal variation of the seismic activity seems to overlap. Multiple seismic events occur within a similar timeframe. An M-T diagram would be helpful for understanding the seismic behavior.

R. Thanks for your suggestion. We have added the M-T diagram in the revised manuscript. The M-T diagram shows that microseismic activities have gradually increased since 2021. This heightened seismic activity is likely indicative of continuous magma movement within the middle-lower crustal reservoirs. 

#Reviewer2

Lü et al. presents a high-resolution 3D shear-wave velocity model beneath the Wudalianchi (WDLC) volcano in northeast China using full-wave ambient noise tomography. They analyzed data from 30 broadband seismic stations deployed around the WDLC volcano to image the crustal structure with improved resolution compared to previous studies. Their key finding is a distinct low-velocity body at depths of ~10-23 km, interpreted as a middle-lower crustal magma reservoir. The data, methods, and results appear robust, and the manuscript is well-written. I have some minor comments.

1. Similar results (low-velocity zones beneath the Songliao Basin and beneath the WDLC volcano) have already been obtained by Lü et al. (2024b). While the authors mentioned in the Introduction section that the resolution is lower than in the present study, there does not seem to be a significant difference overall. Additional explanation should be provided to establish the novelty of this research more clearly.

R. Thanks for your suggestion. The improved resolution model reveals detailed features of the low-velocity zones, providing insights into the geometry of the crustal magma reservoir, a more accurate estimation of the melt fraction, and a clearer understanding of its spatial correlation with seismic activity. Text has been updated in the revised manuscript.

2. Given that WDLC had its last eruption in the 18th century and is considered an active volcano, it would be beneficial to discuss more explicitly the implications for volcanic hazard assessment.

R. Thanks for your suggestion. We assess potential volcanic hazards by integrating the temporal activity patterns of microseismic events with the spatial correlation of low-velocity anomalies. The observed seismic activity, particularly above the middle-lower crustal magma reservoir, may serve as an important precursor to future volcanic unrest. Since 2021, the frequency of microseismic events has shown a notable increase, indicating a progressive intensification of magmatic activity beneath the volcano. This increased seismicity likely points to ongoing magma migration and the pressurization of the middle-lower crustal reservoirs. Therefore, sustained monitoring is essential for detecting precursory signals of potential eruptions, understanding the dynamics of magmatic activity, and mitigating volcanic hazards. We have expanded the discussion in the revised manuscript.

3. As shown in Figures 8 and S5, while the low-velocity zone beneath the Songliao Basin has also been estimated in previous studies, the low-velocity zone directly beneath the WDLC volcano appears unclear in previous study results. Section 4.1 primarily focuses on discussing differences in methodology, but could the authors also mention the contribution of their use of a larger number of seismometers?

R. Thanks for your suggestion. The dense array provides better coverage and higher sensitivity to small-scale velocity anomalies, enabling us to clearly delineate the low-velocity zone directly beneath the WDLC volcano. Text has been updated in the revised manuscript.

4. Could there be further discussion about the time scale and supply volume of magma transport processes between each reservoir in the multi-level magma system (Figure 9)?

R. Thanks for your suggestion. Individual magma reservoir in the upper crust reaches local thermal equilibrium within 10³ years. In contrast, large igneous systems could persist for 10⁵ to 10⁷ years (Cashman et al., 2017). Given that the last eruption of the WDLC volcanoes was ~300 years ago, the upper crustal magma reservoir may be continuously supplied by the observed middle-lower crustal magma reservoir. K–Ar dating results reveal that the Wudalianchi volcano has erupted about 30 times in the past 1 Ma (Wang & Chen, 2005), which preliminarily implies that magma transport to a shallow reservoir may require ~30 ka. In the future, more geochemical and modeling approaches are required to further evaluate magma transport processes in the multi-level magma system. Text has been updated in the revised manuscript.

#Reviewer3

Lû et al. target the shallow-mid crustal magmatic systems beneath the Wudalianchi Volcano in NE China, interesting intraplate volcanism with significant debate on the source origin. To do this, they adopt full-wave ambient noise tomography to seismic data from ~30 stations in the region and derive 3D shear-wave velocity structures from the surface down to 25 km depth. The main finding is a generally low-velocity body beneath the volcano at 10-23 km depth and was interpreted as a mid-crustal magmatic reservoir driving the volcano's activity. The manuscript is well-written, and high-quality figures are used. However, I have a few significant comments regarding the methods and discussion for the authors to consider when preparing their revisions.

L34: what is the context for the statement of “as one of the youngest volcanoes”? Is this in NE China, across China or something else?

R. Corrected.

Section 2 and Figure 3: For the cross-correlation functions, how are the waveforms fitted in the inversion? Are there any specific time windows targeted, or are all correlation functions? At the low-frequency end, 3-7 s and 5-10s, there appear to be multiple phases in the cross-correlation functions (Figures 3a-3b). How are the complex phases handled in the fitting process? In addition, please include waveform fitting plots to show the general quality of their model.

R. Thanks for your suggestion. During the inversion, the cross-correlation coefficient between the synthetics and observed EGFs is required to be at least 0.6, and the observed EGFs are selected with a signal-to-noise ratio threshold of 4. The targeted time windows of Rayleigh waves are calculated with their typical minimum and maximum shear velocity of 2.5 km/s and 5 km/s. We minimize the discrepancies between the synthetics and EGFs with amplitude normalization applied to account for amplitude fluctuations. We have added waveform fitting in the Supporting Information. Text has been updated in the revised manuscript.

Section 3: the authors conducted checkerboard tests and one specific synthetic test to assess the resolution and recovery of their data and methods. These are very useful results, but they tell little information about the vertical smearing effects of their inversion. I suggest adding one or two synthetic tests, similar to the test shown in Figure S3, but putting the low-velocity anomalies at both shallower and deeper depths. The new tests would offer more insights into the robustness of their depth constraints on the mid-crustal reservoir, which is very important information for their model interpretations.  

R. Thanks for your suggestion. We have conducted a new synthetic test to assess the resolution. Figure has been updated in the revised manuscript.

Section 4.1: the comparison plot summarized in Figure 8 is nice. However, it lacks a context of the existing studies. Could the authors add a paragraph to summarize the differences in methods, data used and other significant differences? In addition, the authors should include their own 2024 study for comparison as well. It appears that their 2024 paper uses the same methods but even more stations compared to this study in the region. A detailed justification is required to see significant improvement for a separate study in this paper.

R. Thanks for your suggestion. In comparison with recent shear-wave velocity models from traditional ambient noise tomography (Figs. 8 and S8), such as Shen et al. (2016), Fan et al. (2021), Chen et al. (2021), and Song et al. (2025). The traditional ambient noise tomography is based on the ray-based inversion method and assumes that the Rayleigh wave phase velocity is mainly sensitive to shear wave velocity. Due to wavefront healing (Nolet & Dahlen, 2000), the recovered shear-wave anomalies could be decreased by a factor of two or more (Maguire et al., 2022). Additionally, our deployment of a denser seismic network yields a higher crustal resolution compared to previous regional models (Shen et al., 2016; Fan et al., 2021).  Moreover, compared to the model of Lü et al. (2024b), we retrieve higher-frequency signals that enhance the sensitivity to shallow structure. Text has been updated in the revised manuscript.

Section 4.3: in the schematic model, the authors included the shallow-crust reservoir imaged by a previous study of Li et al. (2016). The authors should consider including a paragraph or two to explain why such a shallow reservoir is missing in their own model. Is this due to the data used or some of the wave complexities containing such information in this study, but they are not considered in the inversion?

R. Thanks for your suggestion. The signal-to-noise ratio of the observed EGFs is relatively lower at short periods due to the wave complexity. Consequently, the shallow reservoir is not well defined due to the sparse ray-path coverage at short periods. Furthermore, with a lateral resolution of ~27 km, our model cannot confidently image small-scale crustal reservoirs with a scale of a few kilometers as reported by Li et al. (2016). Text has been updated in the revised manuscript.