Observations and Modeling of Tsunamis from Genesis to Inundation Using Contemporary Methods by Sean R Santellanes A dissertation accepted and approved in partial fulfillment of the requirements for the degree of Doctor of Philosophy in Earth Sciences Dissertation Committee: Diego Melgar, Chair and Advisor Valerie Sahakian, Core Member Dave Sutherland, Core Member Johnny Ryan, Institutional Representative University of Oregon Fall 2024 c© 2024 Sean R Santellanes This work is openly licensed via CC BY-NC-ND 4.0. 2 DISSERTATION ABSTRACT Sean R Santellanes Doctor of Philosophy in Earth Sciences Fall 2024 Title: Observations and Modeling of Tsunamis from Genesis to Inundation Using Contemporary Methods This dissertation explores broadly the science of tsunamis from genesis to inundation. Tsunami science has been able to produce accurate forecasts of tsunami first arrival and propagation from advancements in deep ocean monitoring and high resolution bathymetry and topography. However, it remains in its infancy compared to its sister disciplne — seismology, owing to the lack and sparse instrumentation in the offshore, deep ocean environment. In addition, the lack of understanding of how turbulence and anthropogenic structures affects inundation leaves the most destructive part of the tsunami life cycle remarkably unconstrained. This dissertation seeks to understand and better constrain several aspects of the tsunami life cycle. In this dissertation, I present an analysis of the tsunami source of the 2022 Hunga Tonga global tsunami, showing from open ocean and coastal sea level stations, the coupled atmospheric conditions that propagated the tsunami globally. Further, I show that the background open ocean tsunami spectrum, vital for understanding tsunami amplification from open ocean to the coast, can experience changes due to atmospherically induced infragravity waves. Next, I present an analysis of the 2020 Sand Point tsunami source. This work shows a novel method for including open ocean instrumentation with geodetic and seismological data to produce finite fault models capable of detecting ”hidden” subduction zone events. Finally, I present an analysis of the differences of inundation between homogeneous and heterogeneous tsunami sources along the Cascadia Subduction Zone. This work highlights the importance of using heterogeneous tsunami sources for understanding how inundation may affect the anthropogenic environment. This work may lay the groundwork for probabilistic- based tsunami evacuation rather than relying on deterministic, homogeneous tsunami sources. This dissertation includes previously published (unpublished) co-authored material. 3 CURRICULUM VITAE NAME OF AUTHOR: Sean R Santellanes GRADUATE AND UNDERGRADUATE SCHOOLS ATTENDED: University of Oregon, Eugene, OR, USA Pennsylvania State University, State College, PA, USA The University of Texas at Austin, Austin, TX, USA DEGREES AWARDED: Doctor of Philosophy, Earth Sciences, 2024, University of Oregon Master of Science, Meteorology and Atmospheric Science, 2020, Pennsylvania State University Bachelor of Science, Physics, 2017, The University of Texas at Austin AREAS OF SPECIAL INTEREST: Tsunamis Earthquakes Data Science PROFESSIONAL EXPERIENCE: Graduate Employee, University of Oregon, 2019-2024 Year-Round Intern, Sandia National Laboratories, 2022-2024 Alumni Volunteer, GeoFORCE Texas, The University of Texas at Austin, 2023-2024 Intern, United States Geological Survey, 2014-2015 GRANTS, AWARDS AND HONORS: PUBLICATIONS: 4 Santellanes, S. R., Young, G. S., Stensrud, D. J., Kumjian, M. R., & Pan, Y. (2021). Environmental conditions associated with horizontal convective rolls, cellular convection, and no organized circulations. Monthly Weather Review, 149(5), 1305-1316. Santellanes, S. R., Ruiz-Angulo, A., & Melgar, D. (2023). Tsunami waveform stacking and complex tsunami forcings from the Hunga-Tonga eruption. Pure and Applied Geophysics, 180(6), 1861-1875. 5 ACKNOWLEDGEMENTS The work done in this dissertation would not have been possible without the help and advice of many people. First, I want to thank my advisor, Diego Melgar, for his endless patience and guidance with the various stuff that came from life and science during the course of this degree. Next, I would like to thank my committee and co-authors Valerie Sahakian, Dave Sutherland, and Dara Goldberg for their help and patience, especially with helping a meteorologist understanding the conventions of seismology. I would like to thank my friends and colleagues that I have made over the course of my degree, specifically, Avigyan Chatterjee, Avery Conner, Kate Scholz, Christina Cauly, Sydney Dybing, Sandia Infrasound Crew, Walker Building 406, and Big Timber Running Club. Additionally, I would like to thank friends and colleagues in the Melgar group, other graduate students, and the Earth Sci front office. And lastly, Sitka for being a good emotional support dog. 6 Dedicated to friends and family who died along the way. . . 7 TABLE OF CONTENTS Chapter Page I. INTRODUCTION . . . . . . . . . . . . . . . . . . . . . . . . 25 II. TSUNAMI WAVEFORM STACKING AND COMPLEX TSUNAMI FORCINGS FROM THE HUNGA-TONGA ERUPTION . . . . . . . . . . . . . . . . . . . . . . . . . . . 29 2.1. Introduction . . . . . . . . . . . . . . . . . . . . . . . . . 31 2.2. Data and Methods . . . . . . . . . . . . . . . . . . . . . . . 37 2.3. Results and Discussion . . . . . . . . . . . . . . . . . . . . . 45 2.3.1. Moveout at 315 m/s . . . . . . . . . . . . . . . . . . . 47 2.3.2. Moveout at 205 m/s . . . . . . . . . . . . . . . . . . . 49 2.3.3. Moveout at 155 m/s . . . . . . . . . . . . . . . . . . . 51 2.3.4. Other effects and observations . . . . . . . . . . . . . . 54 2.3.5. Caveats to the interpretation . . . . . . . . . . . . . . . 54 2.3.6. Hazards implications . . . . . . . . . . . . . . . . . . 55 2.4. Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . 57 2.5. Data and Resources . . . . . . . . . . . . . . . . . . . . . . 58 2.6. Declaration of competing interests . . . . . . . . . . . . . . . . 58 2.7. Acknowledgements . . . . . . . . . . . . . . . . . . . . . . 58 III. CONTEMPORARY MEASUREMENTS OF THE BACKGROUND OPEN OCEAN TSUNAMI SPECTRUM USING THE DEEP-OCEAN ASSESSMENT AND REPORTING OF TSUNAMI (DART) STATIONS . . . . . . . . . . 61 3.1. Introduction . . . . . . . . . . . . . . . . . . . . . . . . . 62 3.2. Data and Methods . . . . . . . . . . . . . . . . . . . . . . . 66 8 Chapter Page 3.2.1. The DART station system and PPSD measurements of the BOOTS . . . . . . . . . . . . . . . 66 3.2.2. Measuring temporal variations in the BOOTS . . . . . . . 67 3.2.3. Measurement of IGWs with DART stations . . . . . . . . 69 3.3. Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . 70 3.3.1. PPSD behavior of the BOOTS . . . . . . . . . . . . . . 70 3.3.2. Spatiotemporal variations in the BOOTS . . . . . . . . . 72 3.3.2.1. Temporal variations in the BOOTS . . . . . . . . 73 3.3.2.2. Spatial variations of power in the BOOTS . . . . . 75 3.4. Infragravity waves in the BOOTS . . . . . . . . . . . . . . . . 78 3.4.1. Spatiotemporal Effects of IGWs . . . . . . . . . . . . . . 78 3.5. Data artifacts in the BOOTS . . . . . . . . . . . . . . . . . . 80 3.5.1. Tsunami signals . . . . . . . . . . . . . . . . . . . . 80 3.5.2. Bifurcation of PPSD in short periods . . . . . . . . . . . 83 3.5.3. Meteorological impacts on the BOOTS . . . . . . . . . . 85 3.6. Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 89 3.6.1. The magnitude of the BOOTS slope . . . . . . . . . . . . 89 3.6.2. Implications of meteorological impacts on the BOOTS . . . . 90 3.6.3. Implications of a spatiotemporally varying BOOTS . . . . . 92 3.7. Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . 93 3.8. Data availability . . . . . . . . . . . . . . . . . . . . . . . 94 3.9. Acknowledgements . . . . . . . . . . . . . . . . . . . . . . 95 IV. AN UNEXPLAINED TSUNAMI: WAS THERE MEGATHRUST SLIP DURING THE 2020 MW7.6 SAND POINT, ALASKA, EARTHQUAKE? . . . . . . . . . . . . . . . . 97 4.1. Introduction . . . . . . . . . . . . . . . . . . . . . . . . . 98 9 Chapter Page 4.2. The incompatibility of a strike-slip source . . . . . . . . . . . . 101 4.2.1. USGS model . . . . . . . . . . . . . . . . . . . . . . 101 4.2.2. Joint model of on-shore and off-shore observations . . . . . . 105 4.3. Deformation requirements of the observed tsunami . . . . . . . . . 107 4.4. Allowing coeval megathrust and strike slip rupture . . . . . . . . 111 4.4.1. Megathrust rupture initiation . . . . . . . . . . . . . . . 111 4.4.2. Full Kinematic Rupture Modeling . . . . . . . . . . . . . 114 4.5. Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 121 4.5.1. Hidden megathrust rupture . . . . . . . . . . . . . . . 121 4.5.2. Potential Submarine Landslide . . . . . . . . . . . . . . 123 4.5.3. Regional context and hazard considerations . . . . . . . . 124 4.6. Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . 129 4.7. Acknowledgements . . . . . . . . . . . . . . . . . . . . . . 130 4.8. Data and code availability . . . . . . . . . . . . . . . . . . . 130 V. COMPARISONS OF TSUNAMI INUNDATION BETWEEN HOMOGENEOUS AND HETEROGENEOUS EARTHQUAKE SOURCES AT SELECT SITES FOR THE CASCADIA SUBDUCTION ZONE. . . . . . . . . . . . . . . . . 133 5.1. Motivation . . . . . . . . . . . . . . . . . . . . . . . . . . 133 5.1.1. Source complexity and coupling in local tsunami hazard estimates . . . . . . . . . . . . . . . . . . . . 137 5.1.2. Are the ASCE 7-22 and DOGAMI ”t-shirt” models still reasonable sources for hazards estimates? . . . . 139 5.2. Data and Methods . . . . . . . . . . . . . . . . . . . . . . . 142 5.2.1. Stochastic modeling of geodetically constrained ruptures . . . 142 5.2.2. Tsunami modeling . . . . . . . . . . . . . . . . . . . 144 5.2.3. Comparisons of hazards between homogeneous and heterogeneous ruptures . . . . . . . . . . . . . . . 146 10 Chapter Page 5.3. Results . . . . . . . . . . . . . . . . . . . . . . . . . . . . 148 5.3.1. Mean tsunami inundation behavior . . . . . . . . . . . . 149 5.3.2. Differences between homogeneous and heterogeneous ruptures. . . . . . . . . . . . . . . . . . 153 5.3.3. Tsunami hazard maps . . . . . . . . . . . . . . . . . . 156 5.4. Discussion . . . . . . . . . . . . . . . . . . . . . . . . . . 157 5.4.1. Differences of homogeneous and heterogeneous models . . . . 157 5.4.2. Moving towards a probabilistic maximum considered tsunami . . . . . . . . . . . . . . . . . . . 161 5.5. Conclusion . . . . . . . . . . . . . . . . . . . . . . . . . . 163 VI. CONCLUSION . . . . . . . . . . . . . . . . . . . . . . . . . 165 REFERENCES CITED . . . . . . . . . . . . . . . . . . . . . . . . 168 11 LIST OF FIGURES Figure Page 1. DART buoys and tide gauges analyzed in this work. DART buoys that were operable with data from only the Tohoku tsunami are shown in neon green. DART buoys that were operable with data for both Tohoku and Tonga are shown in dark orange. DART buoys that were operable with data for only the Tonga tsunami are shown in white. Tide gauges used in this analysis are shown in red. The focal mechanism of the Mw 9.1 Tohoku earthquake is plotted along with the location of the HTHH volcano in white. Atmospheric barometer stations are plotted in teal. . . . . . . . . . . . . . . . 30 2. Schematic diagram showing the 4 potential tsunamigenic source mechanisms from the eruption of HTHH. From left to right, a flank collapse or explosion would have a tsunami propagation speed of ∼220 m s-1. An internal gravity wave would have a propagation speed of ∼150 m s-1. A Lamb wave would have a propagation speed of ∼310 m s-1. Lastly, an acoustic wave would have a propagation speed of ∼330 m s-1. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 35 3. (a) First arrival amplitudes and peak amplitudes in meters recorded at 22 DART buoys from the M9 Tohoku-Oki tsunami and 34 DART buoys from HTHH. (b) Peak amplitudes in meters recorded at 27 tide gauges from the M9 Tohoku-Oki tsunami and 82 tide gauges from HTHH with clear and obvious tsunami signatures. The 403 tide gauges from Carvajal, Sepúlveda, Gubler, and Garreaud (2022) are plotted to compare our dataset with their results. Only first arrival amplitudes are plotted. . . . . . . . . . . . . . . 36 12 Figure Page 4. a.) Record section of the DART data from the HTHH tsunami. Red line denotes the speed of the first arrival seen at the DARTs. The blue line denotes the speed of potential secondary arrivals. b.) Record section of tide gauge data from the HTHH tsunami. The red line denotes the hypothetical Lamb wave arrival at the tide gauge stations. The blue line denotes the second arrival potentially seen the tide gauge stations. The data are shown for as long as 15 s and 1 min data were available. The tsunami amplitude waveforms are scaled individually to better appreciate the arrivals. . . . 39 5. Amplitude stack plots with time delays in hours. The plots show the concentration of amplitudes in the 7 azimuth gates and full azimuth range with available DART data. 305 m/s moveout speeds are discernible in all but one – 135-180. 155 m/s moveout speeds are discernible in 3: 90-135, 135- 180, 225-270. Moveout speeds in typical tsunami velocities are observable in all azimuth gates. The amplitudes are normalized by the number of DARTs in their respective azimuth gates. . . . . . . . . . . . . . . . . . . . . . . . . . . 42 6. Amplitude stack plots with time delays in hours. The plots show the concentration of amplitudes in the 7 azimuth gates and full azimuth range with available DART data. 305 m/s moveout speeds are discernible in all but one – 135-180. 155 m/s moveout speeds are discernible in 3: 90-135, 135- 180, 225-270. Moveout speeds in typical tsunami velocities are observable in all azimuth gates. The amplitudes are normalized by the number of DARTs in their respective azimuth gates. . . . . . . . . . . . . . . . . . . . . . . . . . . 44 7. a.) The amplitude stack plot of normalized DART waveforms combined assuming a velocity of 215 m s-1. b.) The amplitude stack plot of normalized DART waveforms combined assuming a velocity of 305 m s-1. c.) Peak amplitudes at toffset=0 for the averaged, normalized DART stacked waveforms for velocities from 100 - 400 m s-1 for various azimuth ranges. All azimuths contain 35 stations, 0-45 11 stations, 270-315 8 stations, and 315-360 6 stations. d.) The same as c.) but for tide gauge data. There are 114 stations for all azimuths, 0-45 1 station, 270-315 7 stations, and 315-360 14 stations. . . . . . . . . . . . . . . . . . . . . . . 46 13 Figure Page 8. Examples of multi-taper spectral analysis Wigner-Ville spectrograms of equidistant DART buoys. a.) and b.) are about 12 from their respective sources. C.), d.) and e.) are about 95-121 from their respective sources. A.) and c.) show the dispersion of tsunami waves as they move away from the earthquake source. B.), d.) and e.) show how little dispersion there is with the tsunami source from HTHH with the associated Lamb wave front. 21413 is 700 NM ESE of Tokyo, JP. 32412 is 1000 NM W of Lima, Peru. 51425 is 370 NM NW of Apia, Samoa. 42409 is 247 NM S of New Orleans, LA, USA. . . . . . . . . . . . . . . . . . . . . . . . . 49 9. First arrival amplitudes and peak amplitudes in meters recorded at 22 DART buoys from the M9 Tohoku-Oki tsunami and 34 DART buoys from HTHH are converted to pressure amplitudes via equations 2 and 3. We compare these calculated values to terrestrial barometer station peak amplitudes from 7 coastal barometer stations. . . . . . . . . . . . . 52 10. Station up-time graphs that show the time spans of data used for each DART in this study. . . . . . . . . . . . . . . . . . 64 11. Crosses mark the locations of the DART stations used in this study. DARTs 21420 and 46407 are highlighted because they are emblematic (21320) and non-emblematic (46407) to the reference power law of ω−2. BOOTS slope values. . . . . . . . . 67 12. a.) Probabilistic Power Spectral Density (PPSD) plot for DART 21420 — emblematic of behavior that follows the reference power law of ω−2. Amplitudes corresponding to noise levels of 1 cm, 0.1 cm, and 0.01 cm are shown as dashed, black lines. Periods corresponding to 120 s, 240 s, and 800 s are delineated by solid, orange lines. b.) PPSD for DART 46407 — emblematic of behavior that does not follow the reference power law of ω−2. Amplitudes corresponding to noise levels of 1 cm, 0.1 cm, and 0.01 cm are shown as dashed, black lines. Periods corresponding to 120 s, 250 s, and 800 s are delineated by solid, orange lines. . . . . . . 71 13. a.) 21420 time-series of BOOTS slope (blue) and intercept (orange). b.) 46407 time-series of BOOTS slope (blue) and intercept (orange). Period of low-quality data, resulting in no values for either parameter, is shown in red. . . . . . . . . . . . . 72 14 Figure Page 14. a.) Mean slope of the BOOTS across the Pacific basin for each station’s entire temporal coverage. Color of the circle corresponds to the mean value of the BOOTS slope. Circle size corresponds to the standard deviation of the slop. b.) Intercept of the BOOTS across the Pacific basin. Color of the circle corresponds to the mean value of the intercept. Circle size corresponds to the standard deviation of the intercept. . . . 74 15. a.) Median power in dB at a period of 120 s. Top panel is median power during JJA, and bottom panel is median power during DJF. b.) Median power in dB at a period of 250 s. c.) Median power in dB at a period of 800 s. . . . . . . . . . . 76 16. Infragravity wave heights for the time periods of JJA (red circles) and DJF (green circles). DART locations are shown as black dots. . . . . . . . . . . . . . . . . . . . . . . . . . . 79 17. The 2020 Simeonof rupture zone from Crowell and Melgar (2020) is shown in black, the 2021 Chignik rupture zone from the USGS-National Earthquake Information Center (NEIC) finite fault model for the event is shown in dark blue, and the rupture area for the July 15, 2023, Sand Point earthquake from the USGS-NEIC finite fault model is shown in aquamarine (U. S. Geological Survey, 2017). The surface projection of the strike-slip plane associated with the 2020 Sand Point earthquake is delineated by a dashed red line. . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 81 18. DART 46407 PPSD. The 2020 M7.8 Simeonof earthquake is shown by a dash-dot-dot magenta line (UTC Day July 22, 2020). The 2020 M7.6 Sand Point earthquake is shown in blue (UTC Day October 19, 2020) and by a dash-dot crimson line (UTC Day October 20, 2020). The 2021 Chignik earthquake is shown by a dashed green line (UTC Day July 21, 2021). . . . . . . . . . . . . . . . . . . . . . . . . 82 15 Figure Page 19. a.) PPSD for DART 46407. High noise location (red), low noise location (dash-dot-dot blue), meteorological event (dash-dot green), and tsunami (July 22, 2021 Chignik earthquake, dash magenta) PSDs are plotted for their respective occurrences. Arrows point to their respective phenomena in panels b.) and c.). b.) Meteorological event that occurred on November 16, 2020, at 0600 UTC which produced noise above the 10th percentile in the BOOTS. Maximum/Composite radar reflectivity from the 0600 UTC High Resolution Rapid Refresh (HRRR) model is shown. c.) Location of DART 46407 at a high noise location (UTC Day May 17, 2011, red x) and a low noise location (UTC day July 7, 2021, black x). . . . . . . . . . . . . . . . . . . . . . . 84 20. a.) ECMWF wave height model at 0.5◦ resolution for February 10, 2022 at 12 UTC. DART 46407 is shown by a black dot. b.) ECMWF model at 0.5◦ resolution for February 10, 2022 at 12 UTC. Geopotential height at 850 mbar is plotted at 50 m contours. 10 m wind is shown by filled contours. c.) PPSD for DART 46407 for the meteorological event signal during August 13, 2021. The signal is shown by the dashed green line. . . . . . . . . . . . . . . 87 21. a.) PPSD for DART 43413. Meteorological event signal for August 13, 2021 is shown by the dashed green line. b.) GEFS model at 0.5◦ resolution for August 13, 2021 at 12 UTC. 10 m wind is shown to visualize the wind field of the tropical systems. Pressure reduced to mean sea level (MSL) pressure is shown by the black contours. Contours are plotted at 2 mbar intervals. DART 43413 is shown by a white dot. The approximate locations of tropical storm Kevin and Hurricane Linda’s cores are shown by white dots. . . . . . . 88 16 Figure Page 22. The study area, offshore of the Alaskan peninsula. Demarcations for the Sanak, Shumagin, and Semidi segments from Liu et al. (2020) are shown with dashed black lines. The 2020 Simeonof rupture zone is shown in black (Crowell & Melgar, 2020), the 2021 Chignik rupture zone is shown in dark blue (USGS Earthquake Hazards Program, 2017), and the July 2023 Sand Point rupture zone is shown in aquamarine. The W-phase centroid moment tensors (WCMT) for Simeonof, Chignik, and July 2023 Sand Point are shown in the same colors as their rupture zones. The surface projection of the USGS-NEIC finite fault plane for the 2020 Sand Point earthquake is delineated by a dashed red line. The epicenter and WCMT for the 2020 Sand Point earthquake are shown in gold. The King Cove (KING) and Sand Point (SAND) coastal sea level stations are shown as red triangles. DART stations are shown as light blue triangles. GNSS station AC12 (yellow square) is shown to have undergone 10 cm of subsidence during the 2020 Sand Point earthquake. The inset shows the locations of the coastal sea level stations in Hawai’i and the outline of the main study area in blue. . . . . . . . . . . . . . . . . . . . 99 23. a.) Map showing the vertical deformation resulting from the USGS-NEIC finite fault model (model U0). The dashed black line from A-A’ is the surface projection of the causative fault plane, as inferred by the USGS-NEIC. The hypocenter (star) and WCMT are shown in gold. b.) USGS-NEIC finite fault solution. c.) the observed (black) and modeled (red) tsunami waveforms at coastal sea level stations. d.) the observed (black) and modeled (red) tsunami waveforms at DART stations. DART 46403 is shown with slightly muted colors to indicate that it was impacted by Rayleigh wave contamination. . . . . . . . . . . . . . 102 17 Figure Page 24. The results for model S1. a.) Map of study area showing seafloor deformation implied by model S1a constrained to USGS-NEIC magnitude. Line A-A’ shows the surface projection of the strike-slip geometry, with cross section A-A’ below showing the modeled slip distribution along the strike-slip geometry. Model S1a is constrained to a magnitude of Mw 7.6. b.) Map of study area showing seafloor deformation implied by model S1b without any magnitude constraint, resulting in slip distribution equivalent to Mw 8.0. Line A-A’ shows the surface projection of the strike-slip geometry, with cross section A-A’ below showing the modeled slip distribution along the strike-slip geometry. c.) Observed tsunami waveforms (black) with resulting synthetic waveforms from models in (a.) and (b.) shown in red and blue, respectively. . . . . . . . . . . . 108 18 Figure Page 25. The hydrodynamic model (model H0) results. a.) Region map with the 2020 Simeonof and 2021 Chignik rupture zones and associated WCMTs shown in black and blue, respectively (Crowell & Melgar, 2020; U. S. Geological Survey, 2017). The surface projection of the USGS-NEIC finite fault plane for the 2020 Sand Point earthquake is delineated by a dashed red line. The epicenter and WCMT for the 2020 Sand Point earthquake are shown in gold. The King cove (KING) and Sand Point (SAND) sea level stations are shown as red triangles. DART stations are shown as light blue triangles. GNSS station AC12 location shown by yellow square. The red dashed line shows the surface trace for the W-Phase nodal plane used in the USGS-NEIC finite fault model. The top inset shows the locations of the coastal sea level stations in Hawai’i and the outline of the main study area in blue. The bottom inset shows the modeled seafloor deformation (model H0). The black outline denotes an area where a suspected submarine landslide may have occurred, based on the classic dipole sea surface deformation pattern. b.) The tsunami waveforms (observations in black, model synthetics in red) from the coastal sea level stations. Modeled waveforms for Hilo and Kawaehai (KAWA) are shifted by 6.38 min and 2.13 min, respectively, for temporal consistency with the observed data. c.) The tsunami waveforms (observations in black, model synthetics in red) from the DART stations. Gray boxes in a.) and b.) outline the portions of the coastal sea level and DART observations used in the tsunami inversion. . . . . . . 110 26. Potential megathrust rupture nucleation points. a.) Map view of potential nucleation points (orange triangles), including event hypocenter (N0) and nine additional potential nucleation points (N1-N9). Preferred nucleation point, N3, is shown in blue along with the direction and best fitting speed of propagation, 1 km/s. The 2020 Simeonof and 2021 Chignik rupture zones are shown in black and blue, respectively (Crowell & Melgar, 2020; U. S. Geological Survey, 2017) b.) The weighted RMSE for static triggering at nucleation point N0 compared to RMSEs of nucleation points N1-N9 for four different rupture velocities. . . . . . . 113 19 Figure Page 27. Model S2 results. a.) Results of model S2a. Dashed black line delineates the portion of the megathrust considered in the inversion. Magenta contour and shaded area shows the ¿1 m rupture patch of the inversion rupture zone. Dashed line A-A’ shows the surface projection of the strike- slip plane, with slip distribution on the strike-slip plane shown below. The maximum allowable rupture speed is 1.25 km/s for the megathrust and 3.00 km/s for the strike- slip plane. b.) Results of preferred model S2b (nucleation point N3 with maximum allowable rupture velocity of 1 km/s). Dashed black line delineates the subfaults of the megathrust allowed to partake in the inversion. Magenta contour and shaded area shows the >1 m rupture patch of the inversion rupture zone. Nucleation point 3 is shown as a blue inverted triangle. A-A’ is the surface projection of the strike-slip plane from the USGS-NEIC WCMT. c.) The observed tsunami waveforms (black) compared to the modeled tsunami waveforms from models S2a and S2b in red and blue, respectively. . . . . . . . . . . . . . . . . . . . . 115 28. The associated moment tensors (MTs) for the strike-slip segment and megathrust segment for the four kinematic inversions and their composite MTs. The Kagan angle between each composite MT and the USGS-NEIC WCMT are given for each inversion model. The USGS-NEIC MT is shown for visual comparison. . . . . . . . . . . . . . . . . . . . . 117 20 Figure Page 29. Model K3 results. a.) Map showing the geographical distribution of slip along the megathrust as well as the strike-slip geometry used in inversion K3. Black star show the hypocenter for the strike-slip and a blue star shows the nucleation point for the megathrust. Red lines indicate the updip edge of the two fault orientations. b.) Smoothed slip distribution and rupture time contours for the strike-slip segment. Small gray arrows indicate rake direction, scaled by amplitude of slip. Black star shows the hypocentral location. c.) Same as b.) but for the megathrust segment. Note that the rupture time contours start 29.5 s later, as we assume delayed slip to nucleation point N3. d.) The source time function for the published USGS-NEIC finite fault product (black; USGS Earthquake Hazards Program, 2017), the strike-slip segment of model K3 (green), the megathrust portion of model K3 (yellow), and the total source time function of model K3 (magenta) e.) Observed (black) and model K3 synthetic (red) tsunami waveforms. . . . . . . . . . . . . 119 30. The observed tsunami waveforms (black) compared to synthetic tsunami waveforms (red) resulting from the combination of kinematic model K3 and a submarine landslide signal obtained from model H0. . . . . . . . . . . . . . . 125 31. The proposed rupture zone for the 2020 Sand Point megathrust is shown in gold. The 2020 Simeonof rupture zone from Crowell and Melgar (2020) is shown in black, the 2021 Chignik rupture zone from the USGS-NEIC finite fault model for the event is shown in dark blue, and the rupture area for the July 15, 2023, Sand Point earthquake from the USGS-NEIC finite fault model is shown in aquamarine. The surface projection of the strike-slip plane associated with the 2020 Sand Point earthquake is delineated by a dashed red line. The King cove (KING) and Sand Point (SAND) sea level stations are shown in red. DART stations are shown in light blue. The amount of subsidence at GNSS station AC12 (yellow square) for the 2020 Sand Point earthquake is shown to be 10 cm. The inset shows the locations of the coastal sea level stations in Hawai’i. The blue box shows the location of the main study area. . . . . . . . . . 128 21 Figure Page 32. Map of the study area. Insets show the areas of interest that were used for modeling of tsunami inundation. Black arrows point from the inset to the locations of the areas of interest on the main map. CSZ is denoted by the red subduction zone symbols. The XXL1 ”t-shirt” vertical deformation model is shown by the polar plot. Contours are plotted at 5 m intervals. . . . . . . . . . . . . . . . . . . . . . . 140 33. Heterogeneous slip and vertical deformation maps for the four geodetic families used in this study. Mw9.1 source ruptures for a.) 1 cm/yr, b) Gamma, c) Gauss, and d) Li geodetic models. Subduction zone interface geometry comes from Hayes et al. (2018). . . . . . . . . . . . . . . . . . . . . . 142 34. Violin plots of the mean tsunami inundation behavior for each geodetic model and all models. Probability distribution functions shows distribution of the mean of all inundated cells for each scenario (see Figures 37 and 38 for approximate inundated cell extent). Red dashed lines delineate the 25th and 75th percentiles. Solid black lines delineates the 50th percentile. Blue line and arrow show the mean value of all inundated grid cells for the XL1 scenario. Likewise, red line and arrow show the mean value of all inundated grid cells for the L1 scenario. The XYZ line and arrow show the mean value of all inundated grid cells for the M3 scenario. Violin plots are shown for a) Ocean Shores, WA; b) Newport, OR; and c) Crescent City, CA. . . . . . . . 150 35. Maps that show the L1 percentile of inundation for each grid cell for a) Ocean Shores, WA; b) Newport, OR; and c) Crescent City, CA. Black lines denote the official tsunami evacuation lines used for by the states of Washington, Oregon, and California. Streets and geographic information tiles are provided by the USGS. . . . . . . . . . . . . . . . . . . 151 36. Maps that show the XL1 percentile of inundation for each grid cell for a) Ocean Shores, WA; b) Newport, OR; and c) Crescent City, CA. Black lines denote the official tsunami evacuation lines used for by the states of Washington, Oregon, and California. Streets and geographic information tiles are provided by the USGS. . . . . . . . . . . . . . . . . . . 152 22 Figure Page 37. Maps that show the L1 homogeneous where ”wet” is ≥ 30 cm and dry is < 30 cm (a)), average of 200 heterogeneous where ”wet” is ≥ 30 cm and dry is < 30 cm (b)), and difference between ”wet” and ”dry” grid cells for Ocean Shores, WA; Newport, OR; and Crescent City, CA. Black lines denote the official tsunami evacuation lines used for by the states of Washington, Oregon, and California. Streets and geographic information tiles are provided by the USGS. . . . . . . 154 38. Maps that show the XL1 homogeneous where ”wet” is ≥ 30 cm and dry is < 30 cm (a)), average of 200 heterogeneous where ”wet” is ≥ 30 cm and dry is < 30 cm (b)), and difference between ”wet” and ”dry” grid cells of a) and b) for Ocean Shores, WA; Newport, OR; and Crescent City, CA. Black lines denote the official tsunami evacuation lines used for by the states of Washington, Oregon, and California. Streets and geographic information tiles are provided by the USGS. . . . . . . . . . . . . . . . . . . . . . . 155 39. Maps that show the 1-in-50, 1-in-100, and 1-in-2475 year exceedences for the L-like heterogeneous sources for a) Ocean Shores, WA; b) Newport, OR; and c) Crescent City, CA. Black lines denote the official tsunami evacuation lines used for by the states of Washington, Oregon, and California. Streets and geographic information tiles are provided by the USGS. . . . . . . . . . . . . . . . . . . . . . . 158 40. Maps that show the 1-in-50, 1-in-100, and 1-in-2475 year exceedences for the XL-like heterogeneous sources for a) Ocean Shores, WA; b) Newport, OR; and c) Crescent City, CA. Black lines denote the official tsunami evacuation lines used for by the states of Washington, Oregon, and California. Streets and geographic information tiles are provided by the USGS. . . . . . . . . . . . . . . . . . . . . . . 159 23 LIST OF TABLES Table Page 1. The attributes of the inversions considered in this study. SS is strike-slip, MT is megathrust, CSLS is coastal sea level station, SGNSS is static GNSS, HRGNSS is high-rate GNSS, STR is strong motion accelerometer. * denotes unweighted RMSE and †denotes weighted RMSE. . . . . . . . . . . . . . . . . 103 2. Jaccard similarity indices for for Ocean Shores, Newport, and Crescent City of the homogeneous and heterogeneous spatial rasters of ”wet” and ”dry” tsunami inundation. . . . . . . . . 156 24 CHAPTER I INTRODUCTION Tsunamis have been among the deadliest natural hazards of the 21st century. Regardless of country development status, tsunamis have been able to exact high amounts of casualties from > 230,000 from the 2004 Sumatra tsunami to the > 22,000 from the 2011 Tohoku-Oki tsunami. Differing reasons have been ascribed as to why these earthquakes exacted such high tolls. This dissertation focuses on three such reasons: (1) lack of constraint of the open-ocean background state, (2) source mechanism complexity (non-seismic and seismic sources), and (3) varying inundation behavior. This dissertation seeks to test the long-held assumption that the background open ocean tsunami spectrum maintains a constant spectral slope of ω−2 for the entire Pacific basin. The majority of this work was done by myself, with Diego Melgar providing grant funding and minor edits prior to submission to the Natural Hazards and Earth System Sciences journal. This assumption had been applied to spectral analysis and tsunami source inversions for over 20 years. This dissertation finds that the background open ocean tsunami spectrum does not hold constant for the entireity of the Pacific basin; instead, it varies seasonally as function of infragravity waves generated from extra-tropical and tropical cyclones from across the Pacific. This dissertation finds that the North Pacific (from the Alaska Aleutians to off shore of the US Pacific northwest) had significant departures from 25 the ω−2 reference slope due to the region being a hotspot for infragravity wave production. From Santellanes, S. R., Ruiz-Angulo, A., Melgar, D. (2023). Tsunami waveform stacking and complex tsunami forcings from the Hunga-Tonga eruption. Pure and Applied Geophysics, 180(6), 1861-1875. This dissertation then examines the tsunami source complexity of non-sesimic tsunami sources by investigating the tsunami source of the 2022 Hunga Tonga global tsunami. The tsunami from this event traveled faster and further , seemingly able to propagate without obstruction from landmasses. This dissertation finds that the tsunami source of the global tsunami was transient atmospheric source that traveled with the Lamb wave produced by the violent eruption of Hunga Tonga. This dissertation also finds evidence of gravity wave shedding and dispersion that amplified tsunami signals at regional and global distances. It also finds that resonance with bathymetry and topography may have led to increased tsunami amplitudes along the East Pacific. Another examination of tsunami source complexity came by way of the 2020 Sand Point Earthquake. Due to the complexity of this earthquake, a global team was headed by me to investigate the tremor. Dara Goldberg and Will Yeck provided USGS data and techniques for analyzing far-field data. Pablo Koch provided the annealing code that was used as the basis for conjoining water level, seismic, and geodetic data. Brendan Crowell provided expertise of the Alaskan Subduction Zone. Jiun-Ting Ling provided advice on how to invert data from 26 water level data. Diego Melgar was the principle investigator and grant provider. Conceptualization and writing was done by myself. All co-authors provided edits for the manuscript in review, along with two anonymous peer reviewers and Benjamin Brooks, who was the USGS internal reviewer. The earthquake was classified as a right-lateral strike-slip rupture by the United States Geological Survey National Earthquake Information Center (USGS-NEIC). However, it produced a sizeable tsunami that propagated all the way to Hawai’i, and produced hazardous coastal conditions akin to a subduction zone earthquake. In order to investigate the source complexity of this event, this dissertation combines water level data from open ocean instruments and coastal sea level stations with seismic and geodetic data to create a method capable of discerning it. It finds that the most likely tsunami source came from a sizeable slow subduction zone interface rupture that was concealed by fast left-lateral strike-slip rupture preceding it. Finally, this dissertation considers the effects of source complexity on tsunami inundation for the US Pacific northwest from potential Cascadia Subduction Zone ruptures. This work was done by myself with Diego Melgar providing guidence and minor edits to grammar. Simple, homogeneous earthquakes have been used by state and federal agencies to constrain the tsunami hazard of the Cascadia Subduction Zone for ∼ 30 years. These models have also been used by public-private partnerships to establish the land-use and building codes of coastal 27 towns throught the region local to the subduction zone. Homogeneous sources have been shown to underestimate the hazard of tsunami inundation. This dissertation shows that heterogeneous ruptures based off of various geodetic locking models produce tsunami inundation with significant spatial extent differences. In addition, it is shown that the 1-in-2475 (∼ 0.004) exceedence risk of heterogeneous sources also extends further inland than homogeneous sources with the same exceedence risk. 28 CHAPTER II TSUNAMI WAVEFORM STACKING AND COMPLEX TSUNAMI FORCINGS FROM THE HUNGA-TONGA ERUPTION From Santellanes, S. R., Ruiz-Angulo, A., Melgar, D. (2023). Tsunami waveform stacking and complex tsunami forcings from the Hunga-Tonga eruption. Pure and Applied Geophysics, 180(6), 1861-1875. 29 Figure 1. DART buoys and tide gauges analyzed in this work. DART buoys that were operable with data from only the Tohoku tsunami are shown in neon green. DART buoys that were operable with data for both Tohoku and Tonga are shown in dark orange. DART buoys that were operable with data for only the Tonga tsunami are shown in white. Tide gauges used in this analysis are shown in red. The focal mechanism of the Mw 9.1 Tohoku earthquake is plotted along with the location of the HTHH volcano in white. Atmospheric barometer stations are plotted in teal. 30 2.1 Introduction On January 15, 2022 at 04:14:45 UTC, the Hunga Tonga Hunga Ha’apai (HTHH) volcano (Figure 1) erupted violently. The Surtseyan-style eruption likely had a volcano explosivity index (VEI) of 5 (Yuen et al., 2022) and produced an ash cloud that rose into the stratosphere. Teleseismic models from the United States Geological Survey (USGS) ascribed to the eruption an equivalent earthquake magnitude of Ms 5.8. The eruption produced acoustic waves in the audible portion of the spectrum recorded and noticed as far away as Anchorage, AK, USA and The Yukon, Canada at distances over 9,700 km away. Pressure waves from the eruption had an observed maximum value of 7 hPa in New Zealand and had observed values of 2 hPa at various locations around the planet, (Carvajal et al., 2022; Heidarzadeh, Gusman, Ishibe, Sabeti, & Šepić, 2022; Kubota, Saito, & Nishida, 2022; Matoza et al., 2022). The eruption also generated a Pacific wide tsunami along with possible meteotsunamis in the Caribbean Sea, the Gulf of Mexico, and the Mediterranean Sea. Preliminary reports show that in the Pacific wave heights were as small as 0.15 m at Jeju Island, South Korea and as high as 2.5 m impacting Vanuatu (Carvajal et al., 2022). At oceanic basins with almost no connectivity from the source region, the National Oceanic and Atmospheric Administration (NOAA) CO-OPs reported tsunami amplitudes on the order of 0.1 m (Figure 2). The tsunami is meaningful because it had amplitudes significant enough to produce minor damages worldwide (Carvajal et al., 2022). This global nature of the event is reminiscent of the last great volcano-triggered tsunami, the 1883 Krakatoa 31 event (Satake et al., 2020). The paroxysmal eruption there had a VEI of 6 and produced widespread tsunami impacts (Harkrider & Press, 1967; Paris et al., 2014; Press & Harkrider, 1966; Yokoyama, 1987). Unlike HTHH, the Krakatoa tsunami had remarkable amplitudes in the near field (41 m reported at Banten,Pararas- Carayannis (2003)). For that event, barographs recorded pressure arrivals around the world as well. That tsunami seems to have been much larger in the near field than HTHH and resulted in significantly more loss of life. Like HTHH, the Krakatoa tsunami was observed worldwide and analog mareograms of it have been studied before (Harkrider & Press, 1967; Press & Harkrider, 1966). The nature of the physical mechanisms leading to tsunamigenesis in 1883 is debated still. Processes at the volcanic edifice (e.g. caldera collapse, pyroclastic flows, and the explosion) have been invoked to explain it (e.g. Nomanbhoy & Satake, 1995), especially because near-field amplitudes were remarkably high. But, it appears that, like with HTHH, Krakatoa had tsunami arrivals in different ocean basins at times before typical tsunami wave speeds would allow, hinting at complex atmosphere/ocean interactions(e.g. Yokoyama, 1987). The HTHH eruption provides a unique dataset. It is the only volcanic tsunami source to be recorded globally in the modern instrument age of tsunami detection networks. In particular, the increase in deep ocean assessment and reporting of tsunamis (DART) buoys over the past two decades, coupled with a large tide gauge network, has allowed for the tsunami to be characterized on a world wide scale (Figure 1). Here we use these networks to perform an 32 observational investigation into what physical processes can be invoked as the source or sources for the HTHH tsunami waves. A survey of past studies associated with volcanically triggered tsunamis reveals that there are potentially 3 categories of them (Figure 2). The first, and simplest, is the possibility of deformation processes at the volcano such as a caldera collapse or other large submarine mass movements. Were this to be the dominant source we would expect large near-field amplitudes (as in 1883 Krakatau), significant dispersion of the tsunami time series (given a comparatively localized source area, (e.g. Saito, Inazu, Miyoshi, & Hino, 2014), and a propagation speed, or “moveout”, of the water-level disturbances according to vtsun = √ gh, where g is the acceleration due to gravity and h the depth of the ocean. For the pacific basin assuming a mean depth of ∼4000-5000 m this leads to vtsun ∼ 200-220 m s-1. The second, is the possibility of the tsunami being generated by the force on the water column generated by the underwater component of the eruption or the subsequent collapse of initially sub-aerial pyroclastic density currents into the water (Heidarzadeh et al., 2022). Here too, we’d expect a significant dispersive character and a moveout velocity given by vtsun ∼ 200-220 m s-1. For both of these types of tsunami sources we would also expect rapid attenuation with distance given the relatively small geographic footprint of the source compared to that of great earthquakes (e.g. Okal & Synolakis, 2003). Lastly, we have the possibility of pressure disturbances in the troposphere 33 produced by the initial explosion and by the eruption column generating the sea-surface disturbances. Here, it is not redistribution of mass or deformation processes at the volcano which displace the water column and create the tsunami. Rather, it is pressure perturbations in the troposphere which trigger an ‘inverted barometer” response of the sea surface, thus creating tsunami waves. This category of atmospherically-driven tsunami source itself can be quite complex and is typically thought to potentially include three kinds of possible perturbations (e.g. Watada & Kanamori, 2010). These are: acoustic waves, internal gravity waves, and Lamb waves, each of these is thought to produce enough pressure anomalies in the troposphere to induce a tsunami response of the sea surface. Acoustic waves are caused by adiabatic compression and expansion, which propagate with a velocity dependent on the speed of sound in the medium (vair ∼ 330 m/s). Gravity waves are those whose restoring force is gravity. These are common in the ocean and atmosphere (i.e. tsunamis and convection overtopping), and propagate much more slowly than the speed of sound. Lamb waves meanwhile, are interface waves, like Rayleigh waves in seismology. The difference here is that instead of being an elastic wave trapped in the solid portion of the Earth at the boundary with the atmosphere, it is an acoustic interface wave trapped in the troposphere at the boundary with the ocean or solid Earth. They have been observed before with earthquake and volcanic sources which excite the troposphere (Arai et al., 2011; Watada & Kanamori, 2010) and propagate from the source with a characteristic speed slightly below or equal to that of sound, vLamb = 310-340 m s-1 (Nishida, 34 Figure 2. Schematic diagram showing the 4 potential tsunamigenic source mechanisms from the eruption of HTHH. From left to right, a flank collapse or explosion would have a tsunami propagation speed of ∼220 m s-1. An internal gravity wave would have a propagation speed of ∼150 m s-1. A Lamb wave would have a propagation speed of ∼310 m s-1. Lastly, an acoustic wave would have a propagation speed of ∼330 m s-1. Kobayashi, & Fukao, 2014; Revelle & Whitaker, 1996). Lamb waves are long- period traveling atmospheric disturbances, uncoupled from the water, but that can have enough of a pressure perturbation to induce uplift or subsidence of the sea-surface and generate a tsunami. For these types of atmospheric sources we would not anticipate significant dispersion in the tsunami first arrivals (Nishida et al., 2014; Revelle & Whitaker, 1996). Additionally, it is well-known that all of these atmospheric waves propagate with relatively low attenuation (Hedlin, Walker, Drob, & de Groot-Hedlin, 2012). In this paper we aim to understand which of these myriad processes contribute to tsunamigenesis. We use hydrodynamic observations from the deep- and shallow-water worldwide networks to attempt to identify these different generation mechanisms. To do this we carry out a moveout stacking analysis, as 35 Figure 3. (a) First arrival amplitudes and peak amplitudes in meters recorded at 22 DART buoys from the M9 Tohoku-Oki tsunami and 34 DART buoys from HTHH. (b) Peak amplitudes in meters recorded at 27 tide gauges from the M9 Tohoku-Oki tsunami and 82 tide gauges from HTHH with clear and obvious tsunami signatures. The 403 tide gauges from Carvajal et al. (2022) are plotted to compare our dataset with their results. Only first arrival amplitudes are plotted. well as use the spectra of the data at DART buoys to analyze dispersive behavior, and study the amplitude decay characteristics of the event. We will conclude that the first arrival at both deep- and shallow-water stations is almost certainly due to an atmospheric Lamb wave, whose period can be measured at ∼40 mins, with some potential small contributions from acoustic waves. Heidarzadeh et al. (2022) reported a period band of 30-60 min for atmospheric waves, which is consistent with our results. Significant visible later arrivals, especially in shallow water, may be a superposition of internal tropospheric gravity modes and tsunami modes from sources in the water column. As a point of comparison, we rely on the DART and tide gauge data for the previous significant and well-recorded Pacific-basin scale event — the tsunami from the M9 Tohoku-Oki earthquake ((e.g. Satake, Wang, & Atwater, 2003), Figure 3). 36 2.2 Data and Methods We use data from 35 deep water sites (DART buoys) and 114 shallow-water stations (tide gauges) for the bulk of our analysis (Bernard & Titov, 2015). The 15 s and 1 min sampled DART data are available from NOAA and GNS Science for the HTHH tsunami. For the Tohoku-oki tsunami, we also use 15 s and 1 min tide gauge data available from NOAA (DOC/NOAA/NOS/CO-OPS ¿ Center for Operational Oceanographic Products and Services, National Ocean Service, NOAA, U.S. Department of Commerce, 2007). The HTHH tide gauge data is available from the Intergovernmental Oceanographic Commission (IOC) at 1 min intervals (Flanders Marine Institute (VLIZ), Belgium & Intergovernmental Oceanographic Commission-UNESCO, France, 2021). The data were de-tided and band-pass filtered from 5 min to 120 min where the tsunami energy is mostly concentrated (e.g. Rabinovich, 1997). Periods shorter than 5 min contain mostly infragravity waves and periods longer than 120 mins are mostly influenced by the tides. We then plot the tsunami waveforms together, offset by great circle distance from HTHH volcano to produce a record section shown in Figure 4. DART data for HTHH is mainly concentrated in the Pacific basin (32 sites). In addition, we utilize 2 DARTs from the Atlantic basin and 1 DART from the Indian basin. In contrast, all the DART stations for the Tohoku-oki earthquake are only in the Pacific basin since the waveforms were not of meaningful amplitude at other oceanic basins. We show only tide gauges that have a signal-to-noise ratio > 5 and no time gaps (i.e., a clear and obvious tsunami signal) (114 stations out of >800). We calculate the 37 arrival times based on the origin time of the eruption from the USGS event page (2022-01-15 04:14:45 UTC). Heidarzadeh et al. (2022) reported an eruption time of 4:15 UTC based on the analysis of nearby seismic data. In Figure 4, we manually picked tsunami phase arrivals and employed simple least-squares to find moveout speeds that aligned waveform features that seemed coherent from visual inspection. We recognize that this is a somewhat ambiguous or qualitative approach, so, to further explore whether there are coherent arrivals at other propagation or moveout speeds we employ a moveout stacking grid search approach which we now describe. Our moveout stacking process is thus: we first propose a candidate propagation or “moveout” velocity, v. Then, each waveform (deep- and shallow-water) is delayed a certain number of samples by calculating the time delay, tdelay, between the eruption origin time and the expected arrivals given this candidate moveout speed: tdelay = d v + toffset (2.1) where d is the distance between the source and the station in meters. We explored a range of possible velocities in the interval from 100 m s-1 to 400 m s-1. toffset, meanwhile, is another free parameter. It is an arbitrary offset time that allows for arrivals whose origin is not exactly at the assumed event origin time, we allow toffset values in the range of 0 to 6 hours. The goal of this is to explore whether coherent arrivals can be see originating after the eruption origin time. tdelay is then the total delay time which is divided by the sample spacing of the DART or tide gauge recording, dt, to obtain an integer number of samples to delay each waveform 38 Figure 4. a.) Record section of the DART data from the HTHH tsunami. Red line denotes the speed of the first arrival seen at the DARTs. The blue line denotes the speed of potential secondary arrivals. b.) Record section of tide gauge data from the HTHH tsunami. The red line denotes the hypothetical Lamb wave arrival at the tide gauge stations. The blue line denotes the second arrival potentially seen the tide gauge stations. The data are shown for as long as 15 s and 1 min data were available. The tsunami amplitude waveforms are scaled individually to better appreciate the arrivals. 39 by. After delaying, individual waveforms are normalized by their peak value and then all waveforms from either the DART or tide gauge data sets are summed together. This creates a single “stacked” waveform for a particular combination of v and toffset which is then divided by the total number of waveforms in the stack. This results in a maximum theoretical amplitude in any given stack of 1 which would reflect perfect coherence across all waveforms. For interpretation of the results we assume that noticeably high stack amplitudes for a particular combination of parameters this reflects a common physical process as the origin of the “arrivals” or of the “phase” at that particular moveout and offset time. The amplitude of the stacks as a function of moveout speed and offset time can be seen in Figure 5. We are interested as well in whether there are any potential azimuthal patterns so we also separated the waveform stacks into 7 45o azimuth “gates” or “bins” (Figures 5 and 6). To study the spectral features of the waveforms we use the MTSPEC suite of codes from Prieto, Parker, and Vernon Iii (2009) to produce power spectra using the multi-taper method and investigate the dispersive behavior of the HTHH and Tohoku-oki tsunamis. We also produce spectrograms using the Wigner-Ville (WV) approach and apply it to the tsunami waveforms. Multi-taper techniques use prolate-spheroidal tapers and traditional Fourier techniques and have been demonstrated to be minimum variance, unbiased estimators of power which are much superior to the traditional periodogram approach Prieto et al. (2009). 40 Additionally, we extract the peak amplitudes of the first arrivals and of the overall waveforms to evaluate their distance decay characteristics (Figure 3). In addition to the DART buoys and tide gauges, and to understand the nature of the deep-water pressure signals, we analyzed for HTHH, we studied 7 coastal atmosphere barometer stations (locations in Figure 1). This data was bandpass filtered between 12 min and 120 min, to keep comparisons parsimonious with those from the water level data, and to ensure that pressure peaks were not influenced by synoptic scale (L > 1000 km e.g. extra-tropical cyclone) and/or mesoscale (meso-beta and meso-gamma i.e., L < 100 km e.g. thunderstorms or sea breezes) weather phenomena. These features were confirmed to be present during the event via satellite observation from GOES-17 and Himawari-8, and may impact pressure fluctuation calculations if not accounted for in the analysis. We use equation 2 as a starting point to see how Lamb waves may interact with the vertical water column: P ′ = ρCsw (2.2) Where P ′ is the maximum excess pressure available to drive a tsunami response of the water, ρ is the ambient air density, Cs is the air sound speed near the sea surface, and w is the vertical velocity of the water column. w = η T , where η is the maximum displacement of the water column in meters and T is the “rise time”, the period of time water column takes to be displaced in seconds. We vary T 41 Figure 5. Amplitude stack plots with time delays in hours. The plots show the concentration of amplitudes in the 7 azimuth gates and full azimuth range with available DART data. 305 m/s moveout speeds are discernible in all but one – 135-180. 155 m/s moveout speeds are discernible in 3: 90-135, 135-180, 225- 270. Moveout speeds in typical tsunami velocities are observable in all azimuth gates. The amplitudes are normalized by the number of DARTs in their respective azimuth gates. 42 from 0.0005 s to 5 s to systematically test what value is necessary to converted the DART data to atmospheric pressure, the goals of this iterative process is to obtain pressure values similar to that recorded at sub-aerial atmospheric pressure sites. In effect, we are assuming that there should be a match between pressure anomalies recorded at atmospheric pressure sites and the atmospheric pressures needed to drive the responses of the DART buoys. We assume ρ = 1.3 kg m-3 and Cs=330 m s-1. Because the phase speeds of the potential Lamb wave are close to the speed of sound in the case of HTHH, we must employ a correction factor to the equation proposed by Watada (2009) of the form: P ′ = ρoCsw N kCs 1√ 1− (Cp Cs )2 (2.3) Where Cp is the phase speed (305 m s-1), N is the buoyancy frequency (∼ 0.019 rad Hz), and k is the horizontal wavenumber(k = 2π CpTλ ), and Tλ is the period of the Lamb wave (Watada, 2009). In the results section we will show that Tλ can be measured as ∼ 40 min from observation of the DART waveforms tacks. Upon assuming such value for Tλ, we arrive at k∼8.33 x10-6 m-1. In seismology terms, the duration that the Lamb wave deforms the vertical water column, the rise time, would be equivalent to the period of deformation the vertical water column experiences by seismic sea surface deformation, e.g., for Tohoku-Oki’s tsunami this was assumed to be ∼ 30s (Arai et al., 2011; Fujii, Satake, Sakai, Shinohara, & Kanazawa, 2011). However, in this case, we are 43 Figure 6. Amplitude stack plots with time delays in hours. The plots show the concentration of amplitudes in the 7 azimuth gates and full azimuth range with available DART data. 305 m/s moveout speeds are discernible in all but one – 135-180. 155 m/s moveout speeds are discernible in 3: 90-135, 135-180, 225- 270. Moveout speeds in typical tsunami velocities are observable in all azimuth gates. The amplitudes are normalized by the number of DARTs in their respective azimuth gates. 44 estimating the period of deformation due to an atmospheric source acting on the top of the water column. 2.3 Results and Discussion Previous work on Krakatoa and Mt. Pinatubo (Harkrider & Press, 1967; Press & Harkrider, 1966; Watada & Kanamori, 2010) already highlighted the potential complexity of air-sea coupling from an explosive eruption. There it was found that such subaerial sources should produce pressure perturbations in the atmospheric internal gravity, traditional tsunami, and Lamb wave modes. Indeed, Figure 4 suggests from pure visual inspection that there are at least 2 and potentially 3 different tsunami arrivals or phases, with different speeds of propagation. Ostensibly, each is associated with a different generation mechanism (e.g. Figure 2). Figure 5 also shows that across all azimuths there are preferential combinations of tsunami propagation speeds and offset times that produce preferentially high stack amplitudes. In Figure 7 we show the values of these stacks for the tide gauges and DARTs for toffset = 0. This represents a single “slice” through the bottom of each subplot in Figure 5 which captures the behavior if one assumes that all tsunami perturbations originate at the time of the HTHH eruption. We do not find compelling evidence that delayed sources are widespread as we will discuss later on. From Figure 7, at the deep ocean stations (e.g., DARTs) the three potential phases are visible; however, the shallow-water tide gauges only stack somewhat coherently (Figure 7b) for one: 155 m s-1. Following we discuss each potential moveout speed and its significance. 45 Figure 7. a.) The amplitude stack plot of normalized DART waveforms combined assuming a velocity of 215 m s-1. b.) The amplitude stack plot of normalized DART waveforms combined assuming a velocity of 305 m s-1. c.) Peak amplitudes at toffset=0 for the averaged, normalized DART stacked waveforms for velocities from 100 - 400 m s-1 for various azimuth ranges. All azimuths contain 35 stations, 0-45 11 stations, 270-315 8 stations, and 315-360 6 stations. d.) The same as c.) but for tide gauge data. There are 114 stations for all azimuths, 0-45 1 station, 270-315 7 stations, and 315-360 14 stations. 46 2.3.1 Moveout at 315 m/s. In aggregate, the most obvious feature in the data set is a definite tsunami first arrival moving out from the source with a propagation speed of ∼315 m s-1 at DART stations and tide gauges in both the Pacific and Atlantic basins (Figure 4). At the DART buoys it has a far more prominent amplitude than at the tide gauges. It is easily seen in the DART record section as highlighted in red in Figure 4a. The waveform stacks confirm this (Figure 5-7) with a prominent peak in the stack amplitude at these speeds. There is no clear azimuthal variation to the preferential propagation speed and the observed value aligns with what is reported for other observations (Carvajal et al., 2022; Kubota et al., 2022; Matoza et al., 2022) of 300- 307 m s-1 for other related phenomena such as infrasound. This simple observation already likely rules out processes at the volcanic edifice (such as caldera collapse) as the source for the first arrival. Rather that the propagation speed is similar to that of sound in the atmosphere and does not correlate to bathymetry, is consistent with that of a Lamb wave. While Lamb waves from volcanic eruptions have been noted before (e.g. Nishida et al., 2014); this eruption would be the first time that one is observed by the modern tsunami instrument networks. The proposed Lamb wave signature managed to be discerned on all DARTs that were in “event” mode in the Pacific and a few in the Atlantic, which captured it despite not being activated in event mode. Of note as well is that the stack amplitude has azimuthal dependence suggesting perhaps that the Lamb wave generation process has a radiation pattern itself and results in higher amplitudes at certain azimuths. The Lamb wave stack 47 amplitudes are highest between 270-360 and appear to be at a minimum between 90-180 (Figure 5). Additionally, that the stack is so coherent and results in a relatively simple pulse-like shape is indicative of the relative simple propagation behavior – we can measure the period of the Lamb wave at ∼ 40 min from Figure 7. We measure the amplitudes of the first tsunami arrivals as well as on the full waveform (Figure 3). We observe that the amplitudes decay like r−1/2, where r is the distance from the source. This too rules out a localized source area at the volcanic edifice, we would expect much faster attenuation. Note that the HTHH decay is similar to that observed during the Tohoku-oki earthquake. We note as well that while the maximum amplitudes of the Lamb wave are smaller than the maximum amplitude of the overall tsunami waveforms, especially the further the stations are from the source, the pattern of distance decay remains the same. Spectral analysis (Figure 8) is consistent with this view as well. For example, during the tsunami from the Tohoku-Oki earthquake the first arrivals had little to no dispersive behavior, it was not until long distances (>90 degrees) that dispersion became readily apparent in the DARTs. The same is true for HTHH with only marginal amounts of dispersion in the spectra (Figure 8). Revelle and Whitaker (1996) hypothesized, that for Lamb waves, some small amounts of dispersion could be due to differences in meteorology and topography along the Lamb. Finally, we note that while this initial Lamb wave arrival can be seen in some of the tide gauge recordings (Figures 3B,4B) it is very small throughout that 48 Figure 8. Examples of multi-taper spectral analysis Wigner-Ville spectrograms of equidistant DART buoys. a.) and b.) are about 12 from their respective sources. C.), d.) and e.) are about 95-121 from their respective sources. A.) and c.) show the dispersion of tsunami waves as they move away from the earthquake source. B.), d.) and e.) show how little dispersion there is with the tsunami source from HTHH with the associated Lamb wave front. 21413 is 700 NM ESE of Tokyo, JP. 32412 is 1000 NM W of Lima, Peru. 51425 is 370 NM NW of Apia, Samoa. 42409 is 247 NM S of New Orleans, LA, USA. observational data set. It is not the dominant tsunami source in the near-shore environment. 2.3.2 Moveout at 205 m/s. The stack amplitude analysis (Figure 7a) also indicates a peak in the stack amplitudes at a moveout speed of ∼205 m s-1 for the deep-water sites and for all azimuths. Interestingly this speed, however, is not immediately obvious in the record sections (Figure 4). The exact moveout speed varies slightly depending on the azimuth hinting at a process dominated by bathymetry and not by atmospheric 49 propagation. Azimuths that cross deeper ocean depths have faster values around 220 m s-1 and shallower depths having values closer to 205 m s-1. These speeds are consistent with those typically expected in the Pacific basin for regular tsunami propagation. We note as well that moveout speeds in this range do not stack coherently at the tide gauges (Figure 7b) apart from those in the 315-360 azimuth range, which sees an obvious peak at 220 m s-1. There are several peaks in the other azimuth bins that are slightly below 200 m s-1, which may be explained by the shallower bathymetry along the path to those tide gauge stations. This moveout speed is where a “traditional” tsunami generated by processes at the edifice would manifest. That we see a peak in the stacks suggests that there might be indeed some contribution to the recordings from that. However, the distance attenuation pattern proportional to r−1/2 in the peak amplitudes observed in Figure 3, especially for the tide gauges, would indicate that the collapse of the volcanic edifice by itself would be insufficient to account for the basin wide effects observed in Figure 7. Potentially then, some other set of tsunamigenic mechanisms which produce waves that propagate a common tsunami speeds is likely involved. For atmospherically driven events (meteotsuamis) Monserrat, Vilibić, and Rabinovich (2006) detailed 3 potential resonances that may be in play during : Proudman resonance, Greenspan resonance, and shelf resonance. For the DART buoys, Proudman resonance is the more plausible source of the tsunami waves at this moveout speed. The stations are in the open ocean at distances where shelf resonance and edge wave noise should be negligible. To 50 generate a tsunami, this result would imply that the atmospheric source was traveling at a speed that matched that of longwave ocean waves in the open ocean. Of which, the 155 m s-1 and 305 m s-1 velocities would be too slow and too fast respectively —the velocities from 200-220 m s-1 are in the correct range. Typical pressures transients needed for meteotsunami processes are of only a few hPa (Monserrat et al., 2006). The HTHH eruption had an initial pressure jump of ∼7 hPa near the source, and distant atmospheric stations were able to record pressure jumps of ∼ 2 hPa. If we use equations 2 and 3 to estimate what atmospheric pressures are needed to force the DART signals to match atmospheric pressures we obtain the results shown in Figure 9. This suggests that it is possible that the atmospheric pressures needed to generate the amplitudes seen at the DARTs are of a similar order to what has already been identified as necessary in past studies of meteotsunamis (Figure 9). 2.3.3 Moveout at 155 m/s. The third possible moveout velocity where we find potentially coherent stacks is slower and less conspicuous. We see from Figure 4a that there is sometimes a secondary arrival that can be observed at many of the DART stations within ∼20 and between 60-80 and for tide gauges between 70-100 from the source. It propagates with a move-out velocity closer to ∼155 m s-1. Figure 7 shows that for times toffset = 0 this speed does not produce any coherent stacks in the DART data yet appears more strongly in the tide gauge data. Theory and modeling (e.g. Watada & Kanamori, 2010) shows that a more slowly traveling atmospheric 51 Figure 9. First arrival amplitudes and peak amplitudes in meters recorded at 22 DART buoys from the M9 Tohoku-Oki tsunami and 34 DART buoys from HTHH are converted to pressure amplitudes via equations 2 and 3. We compare these calculated values to terrestrial barometer station peak amplitudes from 7 coastal barometer stations. 52 internal gravity mode should exist but that it can often be too close to the tsunami mode itself to be easily visible in hydrodynamic data. If the ∼155 ms-1 mode we observe is indeed a manifestation of the atmospheric internal gravity waves, then it has implications that it is likely generated as a result of the air-sea coupled effects of the Lamb wave. This is speculative, we concede, however note that the ∼155 m/s arrivals identified in the DART gauges in Figure 4, if extended back in time do not connect back to t = 0. This suggests that they are generated with some delay from the original source, this is what motivated adding the toffset term in equation 1, however the stack analysis does not conclusively show delayed sources with propagation speeds of ∼155m s-1. Internal gravity wave generation is ubiquitous in the atmosphere. It is possible that the Lamb wave passage displaced air parcels in a process like frontal passage, generating internal gravity waves that emanate from it. This can occur at those periods which coincide with the tsunami band, and indeed internal gravity waves can have periods of up to ∼17 hours (Ruppert Jr et al., 2022). Using this framework can provides a simple explanation for the internal gravity waves generated far from the volcanic source that we see in Figure 4. However, it doesn’t explain their constant speed of propagation. It may be that the constant speed and period of the Lamb wave and its passage are what lead to the near universal speed of the internal gravity waves when generated around the world. Some of these complexities were discussed by Heidarzadeh et al. (2022), where they showed that the characteristics of the tsunami waves recorded on DARTS and tide gauges 53 were significantly different. The mechanism for how this may be is unknown in the tsunami sciences. Analogs from meteorology may be useful as a starting point. 2.3.4 Other effects and observations. Beyond Proudman resonance and its potential contribution to the deep-water records, the long- duration effects of the meteotsunami at tide gauges is another feature of the HTHH tsunami that warrants an explanation. Without a fully coupled atmospheric- oceanic model it is difficult to ascertain its nature, however it may have been enabled by excitation of Greenspan and shelf resonance. These two resonant effects acting would explain the amplification of tsunami amplitudes in certain areas of the Pacific basin. This longer duration of the event compared to traditional tectonic tsunamis (e.g. Kohler, Bowden, Ampuero, & Shi, 2020) is also apparent in deep water; the codas of the HTHH DART records are significantly longer than for the Tohoku case (Figure 8). This is only true in the Pacific, however, as DART 42409 in the Gulf of Mexico see the Lamb wave passage with no significant tsunami coda being generated. 2.3.5 Caveats to the interpretation. It is clear that the HTHH tsunami is complex and has potentially many processes contributing to its generation. We admit that, as hinted in the previous section, a difficulty in a final interpretation of the event is the lack of a fully coupled atmospheric-oceanic propagation model. Here we have relied on observational evidence to distinguish between different mechanisms. We think the basic interpretation is sufficiently compelling, that the first arrival is a Lamb wave 54 propagating in the atmosphere but that secondary arrivals result from complex atmospheric and oceanic sources perhaps aided by different resonant phenomena. But, ultimately designing a model that replicates the observations will be the final proof of this interpretation. For example, while the waveforms stacks in Figure 1 show slower traveling tsunami phases that do not seem to originate at event origin time, the stacking procedure was inconclusive in ascertaining whether this could indeed be justified. A fully coupled model that replicates the observation is missing. Additionally, we note that the potential source processes at the volcanic edifice are extremely complex. The eruptive plume and its interactions with the ocean and the atmosphere, its time history and the energetics of how it excites both are still a subject of vigorous investigation. Uncertainty in this regard also limits the interpretations we have made here. 2.3.6 Hazards implications. The hazard implications of these complex and simultaneously acting effects are important. Typical basin-scale tsunami warning (Bernard & Titov, 2015) relies on tsunami modeling based on earthquake source parameters and assuming long- waves in the water which have propagation speeds controlled by the bathymetry. For the case of HTHH, there is evidence of a complex set of processes with sources in the water column, such as caldera collapse or the explosion itself contributing only for later arrivals and probably with modest energy while the first tsunami perturbation is due to atmospheric phenomena exclusively. As a result, this first tsunami arrives much earlier than expected at far-field measurement stations and 55 shores. Its propagation speed is faster as it is controlled by the acoustic properties of the troposphere, thus, leading to difficulty in forecasting and communicating the threat of the tsunami to the public. Compounding the problem is that while this first arrival is prominent in deep water sites, its contribution to coastal hazards is indubitably very small. So, forecasting amplitudes and arrival times at sites of interest where people and assets might be at risk is doubly challenging. Future efforts to gauge the hazard of future volcanic tsunamis would benefit from incorporating meteorological and infrasound data. It has been shown in ours and recent works (e.g. Heidarzadeh et al., 2022) that meteorological processes may influence the propagation of air-sea coupled waves, as is the case with infrasound waves (Blom, Dannemann, & Marcillo, 2018). While the number of volcanic tsunamis is sparse, it is possible to model atmospheric profiles under varying conditions with global weather models (e.g., the Global Forecasting System – Finite-Volume Cubed-Sphere Dynamical Core (GFS-FV3)). GFS-FV3 serves as the dynamical core for the Goddard Earth Observing System (GEOS) model, which may serve as a starting point to model varying volcanic sources. It may potentially lead to a probabilistic hazard map of volcanic tsunamis. Overall, the coastal effects of HTHH were small, at least compared to great earthquake tsunamis such as the Tohoku-oki event. However, conditions may still allow for tsunamis and tsunami effects to be hazardous, even at oceanic basins disconnected from the original tsunami source and on a worldwide scale. Events 56 like HTHH are rare “edge cases” but modeling improvements coupled with rapid identification of the source as volcanic in nature could abate the problem 2.4 Conclusion Here we have described in detail the data collected from the worldwide tsunami observing network for the HTHH event including 35 DART buoys and 114 tide gauge stations. The eruption of HTHH produced the largest volcanic source tsunami in the modern age of tsunami observation and it was recorded on a global scale. We have identified three potential phases of the tsunami. The first arrival is almost certainly a Lamb wave phase which travels faster than a traditional tsunami. Secondary phases are most likely a complex amalgam of traditional tsunami and atmospheric effects such as internal gravity waves. Additionally, the tsunami codas are long, we posit that the Pacific basin portion of the tsunami is likely to have been the result of processes like those in meteotsunamis such as Proudman resonance induced by the Lamb wave front, with various coastal areas seeing also effects from Greenspan and shelf resonance. Overall, we believe, that, while processes at the volcanic edifice, such as caldera collapse, likely contribute to tsunamigenesis, the evidence does not support this being responsible the predominant source mechanism. The hazards implications of this are important, the tsunami traveled faster than expected making characterizing it and communicating the threat to the publish very challenging. 57 2.5 Data and Resources DART data was accessed from two different API’s: NOAA and GNS New Zealand. Information concerning the GNS New Zealand GeoNet API and DART data availability can be found at the following: https://github.com/GeoNet/ data-tutorials/blob/main/TILDE/README.md . The DOI for GeoNet’s DART is: https://doi.org/10.21420/8tcz-tv02. NOAA DART data is available from the following: http://www.ndbc.noaa.gov/dart data.php. The data used in this study are archived on GitHub (https://github.com/ssantellanes/Hunga Tonga Water Level) and archived on Zenodo (Santellanes et al., 2022). 2.6 Declaration of competing interests The authors declare no competing interests. 2.7 Acknowledgements We acknowledge the New Zealand GeoNet project and its sponsors EQC, GNS Science, LINZ, NEMA and MBIE for providing data/images used in this study. The work of their engineers and software developers allowed for deep insights into their DART data. Their work provided us with extremely high quality and performance data unparalleled in other areas of the DART network. We also acknowledge Josef Dufek and Thomas Giachetti for providing insights into the processes involved with volcanic eruptions. We also acknowledge the Flanders Marine Institute (VLIZ); Intergovernmental Oceanographic Commission (IOC); Sea level station monitoring facility for the tide gauge data as well as the Center 58 https://github.com/GeoNet/data-tutorials/blob/main/TILDE/README.md https://github.com/GeoNet/data-tutorials/blob/main/TILDE/README.md https://doi.org/10.21420/8tcz-tv02 http://www.ndbc.noaa.gov/dart_data.php https://github.com/ssantellanes/Hunga_Tonga_Water_Level https://github.com/ssantellanes/Hunga_Tonga_Water_Level for Operational Oceanographic Products and Services (CO-OPS/NOAA) for the historical tide gauge data. 59 60 CHAPTER III CONTEMPORARY MEASUREMENTS OF THE BACKGROUND OPEN OCEAN TSUNAMI SPECTRUM USING THE DEEP-OCEAN ASSESSMENT AND REPORTING OF TSUNAMI (DART) STATIONS From Santellanes, S. R. & Melgar, D. In review. Contemporary measurements of the background open ocean tsunami spectrum using the Deep- ocean Assessment and Reporting of Tsunami (DART) stations. Natural Hazards and Earth System Sciences. 61 3.1 Introduction Among the challenges affecting tsunami science today is the need for accurate and expedient source models of tsunamigenic earthquakes. Studies have been been undertaken to constrain the aleatory and epistemic uncertainty of the seismic portion of tsunami source models (i.e. via GNSS (Melgar & Bock, 2015)); however, there has been less focus on constraining these uncertainties for the oceanographic portion (i.e, via ocean bottom instruments (Tsushima, Hino, Fujimoto, Tanioka, & Imamura, 2009)). Rabinovich (1997) described a method to derive tsunami source spectra given by the spectral ratio of the observed coastal tsunami spectrum and the atmospheric-induced wave spectrum convolved with the background open ocean spectrum (BOOTS) and a constant. This method is useful because, if correct, it has the ability to allow separation of path (from bathymetry/topography) and source effects. However, its assertion that the BOOTS slope follows — consistently — a reference power law of ω−2, where ω is angular frequency, is itself a potential source of aleatory and epistemic uncertainty. The basis of this reference power law can be found in the interpretations of the bottom pressure recorder (BPR) data for the open ocean for deployments of a few sites that spanned 1-11 months presented by Kulikov, Rabinovich, Spirin, Poole, and Soloviev (1983) and Filloux, Luther, and Chave (1991). Following those findings, the simple log-linear spectral decay shape for the BOOTS has been traditionally used for tsunami studies focusing on the period band 10 mins to 120 mins. However, there is compelling evidence of a 62 myriad of forcings (Aucan & Ardhuin, 2013; Rawat et al., 2014; Webb, Zhang, & Crawford, 1991) that can potentially be impacting these and shorter periods of the BOOTS and potentially introduce deviations from this simple model (Filloux et al., 1991; Kulikov et al., 1983; Rabinovich, 1997; Rabinovich, Candella, & Thomson, 2013; Webb et al., 1991). The reference power law for the BOOTS has presented a precarious situation by ignoring the short periods of the tsunami spectrum. There is no change from oceanographic forcings; infragravity waves do not appear to influence its shape; and most puzzling of all, atmospheric forcings have been assumed to have no impact despite their purported involvement in the shorter periods by earlier literature (Filloux et al., 1991; Kulikov et al., 1983). The two simplest interpretations of these observations are: (1) The BOOTS characterizes an inherent energy cascade from low frequencies to high frequencies regardless of short period forcings within it, or (2) Previous measurements lack the spatiotemporal longevity to discern the impact of these short period forcings on its shape. Although it is surprising that no deviations have been demonstrated from a reference power law of ω−2. It seems certain, from the work done in this study, that variations in the BOOTS are relatively small at these spatiotemporal scales. Here we will find that, if variations in the BOOTS can be observed, we can only infer these characteristics by using instruments that are not confined to short, modest deployments with modest spatial coverage. Fortunately, instruments, such as the Deep-ocean Assessment and Reporting of Tsunamis (DART) stations, have been designed to handle the rough nature of observation in the open ocean to 63 provide information for tsunami early warning systems and record continuously. We leverage these instruments to measure the long-term behavior of the BOOTS on time scales of 1-10+ years, (Figure 10, Figure 11). Figure 10. Station up-time graphs that show the time spans of data used for each DART in this study. The tsunami community continues to make significant progress towards constraining epistemic uncertainty in source models (Mori et al., 2022). However, a potential complication these efforts may face is the lack of constraint of the open ocean background noise necessary to produce tsunami source models. Constraining this noise is necessary for understanding probabilistic effects of tsunamis as they 64 propagate from the open ocean to coastal locations. As pointed out by Rabinovich (1997), the BOOTS slope and intercept are vital for reconstructing tsunami source spectra, which allows for higher fidelity tsunami source models from the DART stations. They can then be used to diagnose — quickly – tsunami hazards and effects. It is imperative that the BOOTS slope be measured to ensure that its effects on coastal locations can be lead to reductions in epistemic uncertainty. Here, we report that the BOOTS slope varies substantially from the ω−2 reference power law used by previous literature. In fact a simple log-linear power law decays is too simplified a view of the long term behavior of the BOOTS. We find that, while the power law can be suitable for some areas of the Pacific basin and for shorter time spans, it performs poorly in other regions (e.g., the Gulf of Alaska and the Cascadia Subduction Zone). We find that the areas that deviate from the reference power law are in areas where infragravity wave (IGW) production is significant during the time periods of December, January, and February (DJF) or June, July, and August (JJA). We further show that these time periods of IGW production correspond most likely to meteorologicaly significant events (e.g., bombogenesis, extra-tropical cyclones, etc.). Additionally, we observe potential IGW effects on observed tsunami signal power from periods of 120 s - 800 s. We also identify several artifacts in the PPSD, most notably a bimodal behavior, or a bifurcation, in power for periods <80 s, in some DART stations. We conclude that this is most likely due to changes in the DART stations themselves 65 (i.e., upgrades of the bottom pressure recorder to DART II or DART 4G). Finally, we discuss the implications of spatiotemporal varying BOOTS slope. 3.2 Data and Methods 3.2.1 The DART station system and PPSD measurements of the BOOTS. To test the hypothesis of the ω−2 reference power advocated by previous literature, we utilize the DART bottom pressure recorder (BPR) quality- controlled data provided by the National Centers for Environmental Information (NCEI) (Deep-Ocean Assessment and Reporting of Tsunamis (DART(R)), 2005). In real-time operations, DART stations operate using event detection to trigger changes in sampling rate. This variability in sample rates is unsuitable for this study, however, the bottom pressure recorder (BPR) component of the DART system continuosly logs 15 s sampled data. This data is periodically retrieved by NOAA oceanographic vessels and archived. We thus limit our analysis to include the time periods where BPR data with a 15 s sampling rate are available (Figure 10). We do not use days where data are at sampling rates > 15 s. We then ingest day long segments of the time series into the PPSD calculation following the approach detailed by McNamara and Buland (2004). This approach is desirable because it does not require any cleaning or removal of artifacts in the data. No filter or removal of tidal signals is applied to the data — data are quality- controlled. We measure the PPSD from the periods of 30 s to 2 hrs, and then calculate the PPSDs for the 34 DART stations shown in Figure 11. We elect to 66 use period in seconds in order to be parsimonious with the various methods applied by referenced literature and this study. 120°E 180° 120°W 60°W 60°S 0° 60°N 46407 21420 Figure 11. Crosses mark the locations of the DART stations used in this study. DARTs 21420 and 46407 are highlighted because they are emblematic (21320) and non-emblematic (46407) to the reference power law of ω−2. BOOTS slope values. 3.2.2 Measuring temporal variations in the BOOTS. We calculate the spectral slope and intercept of the BOOTS in two week intervals with the multitaper spectral power spectral density (MT-PSD) methods detailed in Prieto (2022). We calculate the BOOTS intercept, as it provides a simple single parameter reference for the amount and type of noise present in the sampled time 67 period (Rabinovich, 1997). We choose two week intervals to capture the apparent seasonal variations of the BOOTS and to minimize noise from daily and weekly oceanographic and meteorologic phenomena. We perform a least squares fit of the BOOTS spectral slope between the periods of 12 min - 100 min. We find that, generally speaking, the PPSDs for the DART stations fit on a spectrum between two archetypal behaviors: BOOTS with mostly linear behavior (Figure 12a) with slope values near what the earlier literature (Filloux et al., 1991; Kulikov et al., 1983; Rabinovich, 1997) described (a shape we refer to as ”convex”) and those that experience a large amount of variation and departure from the log-linear model (what we refer to as ”dromedary”) (Figure 12b), particularly between the periods of 120 s - 800 s. The lower bound of 12 min was selected as it is the period when the dromedary shape of the PPSD (see Figure 12b) begins to taper off, and it is outside of the upper limit of 600 s (10 min) in the IGW band described by Aucan and Ardhuin (2013). We select an upper limit of 100 min to minimize any long period noise effects. As with the PPSD method, we use only data when sampling rates are at 15 s. In order to minimize tsunami effects, we do not calculate the BOOTS slope when a tsunami occurs within a two week measuring period and for the next two week period afterwards. Tsunami coda effects can be discernible for days after tsunamigenesis due to bathymetric and coastline effects both regionally (Melgar & Ruiz-Angulo, 2018) and globally(Kohler et al., 2020). We use the tsunami database provided by the National Geophysical Data Center / World Data Service: NCEI/WDS Global Historical Tsunami Database (n.d.) to filter 68 out tsunamis generated by earthquakes ≥ Mw7.0, landslides, and volcanoes. This procedure ensures that we are measuring background noise. 3.2.3 Measurement of IGWs with DART stations. Rawat et al. (2014) utilize numerical models and in situ measurements of DART data to track IGW events across the Pacific basin. There they show via tracking of IGW burst events that IGW production is highest in the area of the CSZ, producing some of the largest IGW values in the Pacific basin. We use the data preparation method from Aucan and Ardhuin (2013) which assumes that bottom pressure recordings at IGW periods are free of pressure effects from wind-driven waves and thus that the DART BPR data can be used to extract IGW heights at those periods. Following their approach, we consider that we are examining free monochromatic waves of wavenumber k and that we can relate the bottom pressure amplitude pb to the surface amplitude elevation a via a transfer function M , which is a function of water depth D, with the equation: pb = aM = a ρg cosh(kD) , (3.1) where ρ is water density and g is gravity acceleration. We relate wavenumber k to wave frequency f via Laplace’s dispersion relation, (2πf)2 = gktanh(kD). We use Prieto (2022)’s multitaper code to produce power spectral densities Fp(f) of the DART BPR data. The BPR data are subject to the same constraints as to when the BOOTS slope and intercept were calculated in section 2.2. The transfer 69 function M is applied to the power spectral densities to obtain the surface elevation spectral density: E(f) = M2Fp(f). (3.2) We then solve for the significant IGW height, which is defined as the partially integrated spectrum: HIG = 4 √∫ fmax fmin E(f)df. (3.3) We choose the same fmin and fmax as Aucan and Ardhuin (2013), setting them to 8.3 × 10−4 Hz and 1.1 × 10−2 Hz, respectively. As we did not remove tsunami signals prior to calculating the significant IGW heights, we remove them and the 6 following days from the dataset. For each station, we calculate the time averaged HIG for the time periods of JJA and DJF. 3.3 Results 3.3.1 PPSD behavior of the BOOTS. The BOOTS intersects with what is traditionally considered the IGW band in the periods from 60 s to 600 s(Aucan & Ardhuin, 2013; Rawat et al., 2014; Webb et al., 1991). Figure 12a shows the PPSD for DART 21420 — located 480 km southeast of Miyazaki-shi, Japan, Figure 11. This station’s behavior is best emblematic of the behavior expected from the ω−2 reference power law. Its mean slope value is 2.03 ± 0.19, and its mean intercept value is 10−4.48±0.35 (Table S1). We observe from Figure 12a that even 70 for this relatively ”well behaved” site the PPSD is not an exact straight line; it can have non-negligible variation of as much as 5-10 dB in the 120 s - 250 s band, giving the PPSD a convex appearance. Figure 12. a.) Probabilistic Power Spectral Density (PPSD) plot for DART 21420 — emblematic of behavior that follows the reference power law of ω−2. Amplitudes corresponding to noise levels of 1 cm, 0.1 cm, and 0.01 cm are shown as dashed, black lines. Periods corresponding to 120 s, 240 s, and 800 s are delineated by solid, orange lines. b.) PPSD for DART 46407 — emblematic of behavior that does not follow the reference power law of ω−2. Amplitudes corresponding to noise levels of 1 cm, 0.1 cm, and 0.01 cm are shown as dashed, black lines. Periods corresponding to 120 s, 250 s, and 800 s are delineated by solid, orange lines. Figure 12b shows the PPSD for DART 46407 — located 390 km to the west of Coos Bay, OR (Figure 11). Its behavior is dromedary in appearance — with variability as high as ∼15 dB variability and hummocky between the periods 120 s - 800 s — and a significant departure from the ω−2 reference power law, in stark contrast of DART 21420 ( Figure 12a). Its mean slope value is 1.90 ± 0.22, and its mean intercept value is 10−4.25±0.40 (Table S1). Critically, Figure 12b indicates that amplitude spread for periods between 120 s - 800 s is considerable between the 71 10th percentile and the 90th percentile, which may be evidence of external forcing. Lastly, DART 46407 has an enigmatic bifurcation in its PPSD for periods < 80 s. Something that is not observed in the PPSD of DART 21420 (12). DART 46407’s PPSD suggests that the ω−2 power law is a poor fit for modeling the BOOTS between the periods of 120 s - 800 s, coincident with the IGW band’s 60 s - 600 s ((Aucan & Ardhuin, 2013)). This deviation from the expected trend has been noted in previous literature, resulting in them only using the reference power law from periods of only 10 mins - 120 mins (Rabinovich, Candella, & Thomson, 2011; Rabinovich & Eblé, 2015). However, it can be seen in Figure 12 that there are possible effects even in the shorter periods of that range, potentially from IGW power leaking into this band. The BOOTS extends from 2 mins - 120 mins — which includes potential IGW interactions — so it is imperative to constrain what these effects are for shorter periods within the BOOTS. 3.3.2 Spatiotemporal variations in the BOOTS. Figure 13. a.) 21420 time-series of BOOTS slope (blue) and intercept (orange). b.) 46407 time-series of BOOTS slope (blue) and intercept (orange). Period of low-quality data, resulting in no values for either parameter, is shown in red. 72 3.3.2.1 Temporal variations in the BOOTS. We observe from Figure 13a that DART 21420 experiences some amount of seasonal variation over the course of its deployment from December 2019 to October 2021. It achieves maximum slope values > 2.3 during the time period of JJA; meanwhile, it obtains minimum slope values < 1.8 during the time period of DJF. We note that there is a seemingly inverse relationship between the maximum/minimum values for the BOOTS slope and intercept. Lower magnitude intercept values have been attributed to periods of less atmospheric noise (Rabinovich, 1997), meaning periods of calm weather correspond to high BOOTS slope values. We note that the period of minimum BOOTS slope values corresponds to extratropical cyclone season for the north Pacific. However, the length of record for DART 21420 is a modest 2 years long compared to other stations we analyze (Figure 10). Figure 13b shows a much longer 12 year time series showcasing rich seasonal variation of bottom pressure for DART 46407. It obtains maximum slope values in JJA > 2.2 and minimum slope values in DJF <1.4. While partially evident in Figure 13a, Figure 13b demonstrates the seemingly inverse nature of the BOOTS slope and intercept over a much longer period of record — January 2010 - April 2022. As discussed earlier, lower magnitude intercept values correspond with less atmospheric noise and vice versa. Figure 13b indicates that JJA is period of seasonally calm atmospheric noise while DJF is a period of disturbed atmospheric noise. Similar to DART 21420, DART 46407’s seasonality appears to coincide with the extratropical cyclone season in the north Pacific, which multiple sources 73 attribute to producing IGW events (Aucan & Ardhuin, 2013; Rabinovich & Eblé, 2015; Rawat et al., 2014; Webb et al., 1991). Figure 14. a.) Mean slope of the BOOTS across the Pacific basin for each station’s entire temporal coverage. Color of the circle corresponds to the mean value of the BOOTS slope. Circle size corresponds to the standard deviation of the slop. b.) Intercept of the BOOTS across the Pacific basin. Color of the circle corresponds to the mean value of the intercept. Circle size corresponds to the standard deviation of the intercept. We demonstrate in Figure 14 that the other DART stations fall somewhere between the previously demonstrated behaviors. We note that sites with behavior similar those similar to DART 21420 are its neighbors in the southwest Pacific. Similarly, DART 46407’s neighbors in the northeast Pacific exhibit similar behaviors. The exception to this trend is DART 51407 — the DART nearest Hawai’i, Figure 11. This station is notable for the sheltering effect it experiences from the Big Island (Rawat et al., 2014; Webb et al., 1991). It is not affected like other DART stations by meteorological events in the north Pacific, naturally 74 filtering out the noise in the BOOTS. However, we find it is susceptible to meteorologically induced IGW events from the south Pacific and the east Pacific — west of Mexico and Central America. 3.3.2.2 Spatial variations of power in the BOOTS. As previously mentioned, DART stations 21420 and 46407 experience varying degrees of seasonal variations with other stations falling on a spectrum between the two. We show that the stations in the north Pacific experience to varying degrees maximum slope values during JJA and minimum slope values during DJF (Figure 15). The trend is reversed for stations from off the coast of Mexico to off the coast of Peru and in the south-central Pacific. For these reasons, for examining the spatial variations of the BOOTS, we split the time periods into DJF and JJA. We further focus our attention to the periods of 120 s (the start of the tsunami band), 250 s (the middle of the dromedary hump), and 800 s (where the dromedary hump abates and end of the IGW band). We see from Figure 15a that the median power at the period of 120 s is in the -55 - -50 dB range for most of the north Pacific. The southwest Pacific is where background noise is at a global minimum for the JJA and DJF months. The region of