     [file DocuSeg.mc.texi not included]

2020, Ch06a--Amcoh.texi, Creative Commons Attribution-NonCommercial-NoDerivatives 4.0 International License, cc--by--nc--nd.texi, cc--by--nc--nd.4.0

-1

History
*******

-2

Tdr.0.4, 2020-May-28, 2020 Christopher Dimech, CC By-Nc-Nd 4.0

Tdr.0.1, 2019-Nov-29 Friday, 2020 Christopher Dimech, CC By-Nc-Nd 4.0, MPhil to PhD - Transfer Report., Report for research project entitled: "Ambient Geotomography within the Central Mediterranean Sea"., Naiad Informatics

6 Ameliorating Waveform Coherency
*********************************

1

/_/ \_\_.__/\__,_|_|_\___\__,_/__/\___/_||_|

6.1, Abductive Reasoning

Abductive Reasoning is a form of reasoning that is based on Misinterpretation (e.g., dubious inferences) or on Parsimonious Schemes. Misinterpretation highlights the problem that wrong inferences can be made even if the inference is based on experimental data.  Parsimony is based on simplicity.  Occam Razor, for instance, adopts the statement that when one is confronted by multiple hypotheses about the same observation, the hypotheses with the fewest assumptions is selected -  a strategy that is also in accordance with English Cosmologist Stephen William Hawking, German Polymath Gottfried Wilhelm Leibniz and English Polymath Isaac Newton Hiroshi & Gibbs (1997)Hiroshi & Gibbs (1997): What is Occam's Razor?.

The Correlation Function as the basic principle by which weak wavetrains in coherent states are identified and extracted from within chaotic wave fields has shown great promise, growing spectacularly in the past three decades Duvall, Jefferies et al. (1993)Duvall, Jefferies et al. (1993): Time-distance helioseismology..  Extremely detailed observational studies based on long, almost uninterrupted recordings of oscillation amplitudes were required to confirm the value of this statistical signal identification tool [tbc06--PractBenf--ClSeism].  This development has come in large part from the meticulous work done on the subject of Wavefield Pre Processing, Cross Correlation Workflows, and the heavy reliance on Stacking, in an effort to substantially improve measurements of Wave Travel Times, estimates of the Green Function of the Seismic Wave Equation, and subsequent Passive Imaging of the Subsurface [tbc06--ProcessingPl--GeoHlSeism].

But as Computational Science and Engineering (CSE) has revolutionised scientific investigations, it is still a relatively new discipline in many Academic Departments.  This is because many professors do not value computational considerations or do not understand the full benefits of computational tools for the improvement on geophysical quantities.  Canadian Mathematician Daniel Lemire and Russian Mathematical Physicist Ludvig Dmitrievich Faddeev have given a similar analysis of the pressure of old authorities on physical thinking and terminology.  Lemire in particular points out that because most professors spent little time programming, good practices such as code reviews and unit testing are not enforced Lemire (2012)Lemire (2012): On the quality of academic software..  Furthermore, most academic software is written by inexperienced students who only write software for their own use and have no idea on how to come up with good designs Cook (2012)Cook (2012): Writing software for someone else.. And although maintenance and documentation are discussed in some courses, they are rarely put into practice.

Faddeev (2000)Faddeev (2000): Modern mathematical physics: What it should be. gives two revealing examples of the extreme attitude and influence of well regarded and established physicists and mathematicians: 1) the criticism by Austrian Theoretical Physicist Wolfgang Ernst Pauli towards Guage Particles for the Quanta of Yang-Mills Field Theory; and 2) the censorship imposed by Soviet Azerbaijani Physicist Lev Davidovich Landau regarding the renormalisations used to treat infinities arising in Quantum Field Theory Yang & Mills (1954)Yang & Mills (1954): Conservation of isotopic spin and isotopic gauge invariance..  Today, Yang Mills Field Theory (also known as Non Abelian Gauge Theory) forms the basis of the Standard Model of Particle Physics, unifying Electromagnetic Forces and Weak Forces; and in the Theory of the Strong Nuclear Force (i.e., in Quantum Chromodynamics) Jaffe & Witten (2009)Jaffe & Witten (2009): Quantum Yang Mills Theory..  See also the complementary comments in Faddeev (2005)Faddeev (2005): Mass in Quantum Yang Mills Theory..  And Renormalization has been shown to be required in one of the pillars of contemporary physics - namely in Statistical Field Theory Weinberg (2013)Weinberg (2013): The quantum theory of fields..

With the novel geophysical pipelines associated with passive imaging NwIdea--PassiveImg, Seismology and Geophysics could change drastically in the next decades, with Geomathematics and Geoinformatics being the driving forces behing new developments.  For this to happen many mathematicians have to become enthusiastic about observational concepts.  Geoscientists, on the other hand, have to become aware of new techniques and perspectives, encourage the development and integration of sophisticated mathematical intuitions and computational frontiers, rather than constraining them to merely function as ancillary services of less value.

For the case of microseismicity, the computation of the Cross Correlation computed as a Time Domain Sliding Inner Product between Station Seismograms has proved to be effective in extracting coherent microseismic surface waves, that are normally hidden within extremely noisy ground motion seismograms.  Notwithstanding the argument above, the Cross Correlation is the most parsimonious description of similarity for ambient surface wavefield retrieval. Schimmel & Paulssen (1997)Schimmel & Paulssen (1997): Noise reduction and detection of weak, coherent signals through phase weighted stacks., in particular, argue how the many inadequacies of standard similarity measures can be resolved by using a non linear stacking technique known as Phase Weighted Stacking, where oscillation phases are used instead of amplitudes.

In this chapter I show that for extremely subtle microseismic signal coherencies, the parsimonious methods that were previously regarded as adequate are no longer acceptable.  And because some may argue that the additional computational burdens are not strictly needed, I will start by introspecting how extremely precise signal correlation estimation methods in local helioseismology have had such a profound effect on fundamental physics, that a modification to the Standard Model of Particle Physics became imperative.

Above all, the rich scientific results obtained from the ongoing microseismic coherence work in the Central Mediterranean Region - using Seismic Station Data in Sicily [South of the Italian Peninsula], the largest island in the Mediterranean Sea - is presented and briefly discussed.  Specifically, two sophisticated Microseismic Coherency Measures are employed 1) Time Frequency Phase Weighted Stacking based on the S Transform; and 2) Time Scale Phase Weighted Stacking using Frames of Continuous Wavelets.  The focus, placed on the Multi-Scale Multi-Resolution Techniques is crucial, because as shall be observed and recognised in the next chapter, the Algorithmic Efficiency of the two Time Frequency Coherency Methods can be improved substantially compared to the conventionally used signal processing methods, for instance, those implemented in ObsPy and MSNoise.  The development of ObsPy is described in Krischer et al. (2015)Krischer et al. (2015): ObsPy: A bridge for seismology into the scientific python ecosystem..

MSNoise suffers from the same problems encountered in ObyPy, particularly, the reliance on ObsPy, Conda and MiniConda; see MSNoise 1.5 Release Notes.  The software also relies on MySQL or SqLite for storing storing configuration bits, waveform metadata, and jobs - a requirements that is considered to be an overkill (i.e., excessive and unreasonable).  However, the most serious problem is the presumption by Lecocq, Caudron & Brenguier (2014)Lecocq, Caudron & Brenguier (2014): MSNoise, a python package for monitoring seismic velocity changes using ambient seismic noise. that MSNoise allows only non-commercial educational and research usage but prohibits Commercial Services without signed agreement with the authors - a view that is inconsistent with the European Union Public License Version 1.1 (EUPL v1.1), under which MSNoise is licensed.  In the next chapter, arguments are given to support the claim of inconsistency.

studies and often sufficient for preprocessing.  Compiled programming language solutions using the central processing unit (CPU) or the graphics processing unit (GPU) are however necessary to apply the most demanding operations to large data volumes in a reasonable amount of time.

6.2, Parsimony over Complexity

Test

2a, The Claerbout Conjecture

Conceptually, computing the Interstation Impulse Response of microseismiocally excited waves is represented by the Green Function, which may be extracted by computing the Time Derivative of the Interstation Cross Correlation -  a simple measure of waveform similarity that permits signal identification Fichtner & Tsai (2019)Fichtner & Tsai (2019): Theoretical foundations of noise interferometry..  The use of Cross Correlation was first conjectured by American Geophysicist Jon F. Claerbout after showing that the reflectivity response of a one dimensional (1-D) layered structure from a surface source and a surface receiver is one side of the autocorrelation of the seismogram recorded at the same surface receiver but transmitted from a source at depth Claerbout (1968)Claerbout (1968): Synthesis of a layered medium from its acoustic transmission response..  The conjecture is appealing because the Travel Time and the Wavelet Shape obtained, are dependent on the subsurface properties of the intervening medium.  Thus, instead of isolating actual excitation sources as is done in Earthquake Seismology, a virtual pair acting as source and receiver are used to map subsurface inhomogeneities in ways similar to reflection seismology, routinely used in the exploration of oil and gas reservoirs.

6.3, Correlograms in Helioseismology

The first application and confirmation of the Claerbout Correlation Conjecture was in Helioseismology, specifically, the introduction of Time Distance Duval, Jefferies et al. (1993)Duval, Jefferies et al. (1993): Time-distance helioseismology..  Discussions on the Helioseismic results follow.

3a, Sunspots and their Interior Structure and Dynamics

Sunspot Imaging

3b, Solar Magnetic Cycle and the Interface Dynamo

Through satellite and computer technology, Astrophysics has developed most rapidly during recent years, making it possible to study the different physical processes taking place in stars.  The discovery of propagating sound waves on the Surface of the Sun Leighton, Noyes & Simon (1962)Leighton, Noyes & Simon (1962): Velocity Fields in the solar atmosphere. I. Preliminary Report. (see also Howard (1967)Howard (1967): Velocity Fields in the solar atmosphere.), and their understanding (i.e., the physical and mathematical explanation) as described by Ulrich (1970)Ulrich (1970): The five minute oscillations on the solar surface. (and independently by Leibacher & Stein (1971)Leibacher & Stein (1971): A new description of the solar five-minute oscillation.) has led to the development of Helioseismology.  The common theme is the Virtual Source Correlation Technique, that provided the first important insight into the Subsurface Structures (i.e., the shallow model) and the Dynamic Nature (e.g., Plasma Downflows) of Sunspots.  To put things in perspective, the Virtual Source Correlation Technique follows the same strategy for mapping the Internal Structure of the Earth using Surface Wave Microseismicity driven by Atmosphere Ocean Seafloor Coupling.

Before focusing on the Correlation Techniques for Earth Microseismicity, an exposition of the computational power of the technique is presented. Besides providing an unprecedented tool into imaging the Solar Interior, the results obtained by Correlation Based Helioseismology have led to fundamental changes in the Standard Model of Particle Physics.  First, results about the structure of sunspots is presented.  Then, the Neutrino Oscillation Experiments are is discussed in some detail.  to convince you on the validity and power of employing novel Correlation Based Techniques to the extraction of Surface Wave Wavetrains for superior Microseismic Subsurface Mapping.

The appearance and disappearance of Sunspots - dark, irregularly shaped areas - are the most notable occurrences on the Surface of the Sun.  Only with the invention of telescopes (e.g., by Italian Astronomer Galileo Galilei [1609, 1610]) was it possible for sunspots to be systematically observed Galilei, Scheiner et al. (2010)Galilei, Scheiner et al. (2010): On sunspots..  Regular observation using instrumental solar data (e.g., using optical telescopes) have mostly been made by professional astronomers with doubtless scientific thoroughness, that can be extended back to the early 17th Century (e.g., [1608] Hans Lippershey; [1611] Johannes Kepler; [1686] Christopher Huygens) King (2011)King (2011): The history of the telescope..  A good exposition is provided by Choudhuri (2015)Choudhuri (2015): Nature's Third Cycle: A story of sunspots., an Indian Astrophysicist familiar with Astrophysical Magneto-Hydro-Dynamics (Mhd.)  in the context of the Solar Magnetic Cycle.

Helioseismic analysis of sunspots have shown a combination of plasma inflows and downflows close to the solar surface, where strong magnetic fields block the flow of heat below the surface.  Thermal blocking making sunspots much cooler than their surroundongs at shallow depths. Sunspots are also associated with energetic Solar Flares and Coronal Mass Ejections.  Sunspot Numbers and their Equatorward Drift of Activity (propagation) are correlated with the Solar Magnetic Field generated by the Solar Magnetic Flux Transport Dynamo Mechanism that follows Magneto-Hydro-Dynamic Laws (i.e., Mhd. Magnetic Laws as well as the laws of Fluid Dynamics).

Sunspots follow an 11 Year Activity Cycle, with their number and location varying as the cycle progresses towards Solar Activity Maximum -  where the number of sunspots peaks.  For identification purposes, the Solar Sunspot Activity Cycle can be synthesised into two temporal periods, each centered around Solar Maximum and Solar Minimum. Sunspots begin to appear around Heliographic Latitudes 45^\circ to 30^\circ during the beginning of the Solar Maximum Period, drifting towards the Solar Equator as the cycle progresses.  The sunspot drift continues and upon reaching the beginning of the Solar Minimum Period, sunspots start to appear around Heliographic Latitudes 15^\circ, continuing their drift equatorward igc06--03--SunspotCycle.  Helioseismic Imaging Techniques lend support to the Interface Dynamo Models, in contrast of earlier Convection Zone Dynamo Models where the Magnetic Dynamo Process was assumed to occur throughout the entire convection zone.

The Solar Dynamo Mechanism defines the behaviour through which ioniosed gaseous plasma states that rotate, convect, and electrically conduct maintain a magnetic field over astronomical time scales. (Plasma is a Magneto Fluid State of Matter defined by unbound positive and negative particles.)  For the Sun, the dynamo depends on the internal shearing of magnetic fields by differential solar rotation that is concentrated within a thin layer - the Tachocline - located at the base of the convective zone, where Helioseismic Inversions (i.e., subsurface estimation techniques that are fundamentally different from the classic scientific observation methods in physics) show strong velocity gradients.  However, Tachocline Shearing is not a crucial assumption in Dynamo Like Models, and can be used for Low Mass Fully Convective Stars where Magnetic Field Generation occurs fully in the convective zone.

A simple mechanism for a Fully Convectice Dynamo occurs when there is enough time for the star to experience the Omega Alpha Effect (\omega\; \alpha Mechanism) - a magnetic flux generator induced by Differential Rotation and Magnetic Flux Twisting that is also characterised by the transformation of a Meridional Magnetic Field (North South Field) into an Azimuthal Magnetic Field (West East Field) and vice versa Love (1999)Love (1999): Reversals and excursions of the geodynamo..  Still, although Magnetic Dynamo Models use Magnetic Induction Equations, the generation of Stellar Magnetic Fields is subtle and intricate because of their dependency on the specific form of plasma flows as well as on their magnitude.  For instance, to convert Kinetic Energy into Magnetic Energy, Kinematic Dynamos require complex turbulent dynamics (chaotic property changes) entailing three dimensional helical plasma flows, irrespective of their magnitude. ([Kinematic Dynamos] are primarily specified by a flow field.  In contrast to Laminar Flow, turbulence prohibits the plasma from travels in layers.  Consequently, mixing of plasma is heavy.)

3c, Solar Neutrino Deficiency Experiments

The most notable scientific breakthrough made by Helioseismology was in solving the Reduced Electron Neutrino Flux reported by the Homestake Gold Mine Experiments Cleveland et al. (1998)Cleveland et al. (1998): Measurement of the solar electron neutrino flux with the Homestake Chlorine Detector., with further confirmations from other Neutrino Measurement Experiments; e.g., by the Super-Kamiokande Experiments Fukuda et al. (1998b)Fukuda et al. (1998b): Measurement of a small atmospheric \nu_\mu / \nu_e ratio. ([The Homestake Experiments] used a tank of PerChloroEthylene to detect radioactive Argon-37 [^37\rm Ar]; following the interaction between Electron Neutrinos [\nu_e] and Chlorine-37 [^37\rm Cl].) The experiments were led by American Chemist Raymond Davis and American Astrophysicist John Norris Bahcall.)

Super-Kamiokande was a Water Cherenkov Neutrino Detector (of recoil electrons from ^8B Solar Neutrinos), located near Hida City (Kamioka Area, North Gifu Prefecture, Central of Honshu Island) at a depth of 1km in the Mozumi Mine under Mt. Ikeno-Yama. ([Mozumi Mine, Cadmium Poisoning] Mozumi Mine is notorious for the Cadmium Mass Poisoning that occurred from 1910 to 1945, when mining operators released cadmium from a processing plant into the local water.  Although the Mozumi Mine Cadmium Poisonings constitutes one of the Four Big Pollution Diseases of Japan, the mines are still in operation and cadmium pollution levels remain high.) Consisting of two concentric cylindrical regions holding Ultrapure Water into which Photo Multiplier Tube Sensors were immersed, Super-K detected Cherenkov Radiation - an optical shockwave emitted when charged leptons particles (e.g., muonons) pass through a dielectric media at speeds exceeding the phase velocity of light in the medium (e.g, the blue glow of an underwater nuclear reactor).

The Inner Detector of Super-K had inward facing Photo Sensors, surrounded by an Outer Detector with outward facing Light Sensors [20cm Diameter] used to distinguish 1) \underline\hboxCosmic Ray Muon Events (muons are electrically charged at a value of -e and considered background sources to the experiments) from 2) \rm \underlineNeutrino\; Events (neutrinos being electrically neutral particles).  Observing neutrinos from all directions, Super-Kamiokande determined that the number of upward going Muon Neutrino Events (which travel long distances inside the earth before being detected) were only half the number of down going events.

3d, Neutrino Oscillations and the Revision of
Particle Physics

The discrepancy in the neutrino counting results created the Solar Neutrino Problem.  Initially, the assumption was that the Standard Model of Particle Physics was correct, whilst the Standard Solar Model was considered problematic.  Adopting this belief, the proposed changes necessary to accomodate the Electron Neutrino Flux Descrepancy focused on reducing the core temperature of the sun to 13 \times 10^6K by modifying 1) the Equation of State; 2) Nuclear Cross Sections; 3) Chemical Composition; 4) Element Abundances; and 5) Radiative Element Opacities (opacities are solar parameters that describe the coupling between radiation and matter in the hot dense interior of the sun).  Yet, Helioseismic Inverse Computations consistently showed that although a depression in the sound speed profile near the solar core existed, the central temperature was constrained to T_c = 15.7 \pm 0.3 \times 10^6K, ruling out the lower core temperature adjustment that would have accounted for the Electron Neutrino Flux Descrepancy [igc06--SolrDelc--1999Dziembowski].  The doubts about the suitability of the Standard Solar Model soon vanished.

Motivated by the ensuing experiments in Canada (the Sadbury Neutrino Observatory Experiment [Sno]), in Italy (Gran Sasso Boron Experiment [Borexino]), and Russia (Soviet-American Gallium Experiment [Sage]), it became evident that something was happening to the neutrinos on their way to the detectors.  The puzzle was solved thanks to the Three Neutrino Lepton Flavour Oscillation Mechanism, whereby Neutrinos undergo a quantum mechanical fluctuation that is determined by the superposition of neutrino wave packet eigenstates.  This mixing allows a neutrino to change state, from one flavour to another, and in so doing, changes in size and mass. ([Neutrino Oscillations] Two examples of neutrino transitions are an Electron Neutrino becoming a Muon Neutrino (\nu_e \rightarrow \nu_\mu ) or a Tau Neutrino (\nu_e \rightarrow \nu_\tau).)

The breakthrough came with the Sudbury Neutrino Observatory [Sno] Ultra Pure Heavy Water Experiments.  Sensitive to all neutrino types - not just to electron neutrinos - Sno Experiments determined there were no more missing neutrinos from the Sun.  Furthermore, the transitions depended on both the distance and time of neutrino propagation between the source and detector.  The two observations clearly confirmed that the neutrino count discrepancies were not caused by problems with the Physical Models of the Sun, but with Elementary Particles Physics.  See confirmation using measurements of the ^8B Solar Neutrinos Ahmad et al. (2001a)Ahmad et al. (2001a): Measurement of charged current interactions produced by 8B solar neutrinos at the Sudbury Neutrino Observatory..

A review of the published results for Sno is provided by Bellerive et al. (2016)Bellerive et al. (2016): The Sudbury Neutrino Observatory. from experimental data collected from 1999 to 2006.  The background includes details of the Three Physics Phases as Sadbury together with the combined analysis of all three phases.  The discussion of Neutrino Oscillations centers on the mass differences, the three Mixing Angles (\theta_12, \theta_13, \theta_23), and the Charge Parity Violation Phase Angle \delta_\rm cp.

The Three Flavour Model of Neutrino Oscillation (f \in \e, \mu,\tau\) relies on a relation associating Flavour Eigenstates \nu_f with their Mass Eigenstates \nu_m using the Pontecorvo Maki Nakagawa Sakata (Pmns) Matrix U_e\mu\tau.

The \\theta_12, \theta_13, \theta_23, \delta_\rm cp\ parameters set the amplitude of Neutrino Flavor Oscillation, whilst the frequency is determined by the differences in the squares of the neutrino masses \Delta m^2_ji = m^2_j - m^2_i.

Bellerive et al. (2016)Bellerive et al. (2016): The Sudbury Neutrino Observatory. also document other physics studies based on Sno Neutrino Expermimental Data to test solar models, neutrino properties, and to search for Extra Terrestrial Neutrinos Sources.  Elaborations on Electron Neutrinos from Charged Current Reaction [Cc] on Deuterium [D_2O] are discussed in Ahmad et al. (2001b)Ahmad et al. (2001b): Measurement of the rate of ve + d -> p + p + e^- interactions produced by 8B solar neutrinos at the Sudbury Neutrino Observatory..  Clearly, the signicantly small counts for Charged Current Reactions (\nu_e + d \rightarrow p + p + e^-) where subatomic particles interact by means of the Weak Force carrying Electrically Charged Bosons W^\pm, provided the quantitative evidence that the observed solar neutrinos were not all Electron Neutrinos \nu. ([Charged Boson Interaction] Weak Force interactions by Chared Bosons involves a transfer of electric charge and weak interaction charge.)

6.4, The Wapenaar-Green Equivalence Principle

The proof of Claerbout Wavefield Coherence Conjecture Claerbout (1968)Claerbout (1968): Synthesis of a layered medium from its acoustic transmission response., extended to 3-D arbitrary inhomogeneous media was provided in Wapenaar (2003)Wapenaar (2003): Short Note: Synthesis of an inhomogeneous medium from its acoustic transmission response., long after its routine practical use in Helioseismology Shapiro, Campillo, Stehly & Ritzwoller (2005)Shapiro, Campillo, Stehly & Ritzwoller (2005): High-resolution surface-wave tomography from ambient seismic noise..

6.5, Seismogram Cross Correlation

The Parsimonious Cross Correlation uses the Occam Razor Heuristic by running a sum of products to determine the similarity between two seismograms (the wavefields, w_a and w_b) recorded at two distinct spatial locations (a and b).  The Parsimonious Approach has been very promising as it allowed fast computations through the application of the Cross Correlation Theorem - an equivalent computation that is described in terms of the Fourier Transforms of the oscillation signals w_a and w_b.

The overbar [overline] of \overline\widehatg(\omega) refers to the Complex Conjugation of the Fourier Coefficients of g(t). See Naughton & Roux et al. (2016)Naughton & Roux et al. (2016): Ambient noise correlations on a mobile, deformable array.; Bensen, Ritzwoller et al. (2007)Bensen, Ritzwoller et al. (2007): Processing seismic ambient noise data to obtain reliable broad-band surface wave dispersion measurements.; Villasenor et al. (2007)Villasenor et al. (2007): Ambient noise surface wave tomography of the Iberian Peninsula: Implications for shallow seismic structure.; Roux & Kuperman et al. (2004)Roux & Kuperman et al. (2004): Extracting coherent wave fronts from acoustic ambient noise in the ocean..

For arbitrary functions g(t) and h(t), the equivalent computation is achieved by taking the Complex Cross Spectrum \widehatR_l(\omega) - the point wise product of spectra after complex conjugating one of the Fourier Transforms - followed by an inverse transform.

Using the Functional Suffix Notation \displaystyle\vphantomu \over v\Gamma_gh(t) to represent the Cross Correlation or Cross-Covariance between two Time Series in the Time Domain, the Cross Spectrum in the Frequency Domain \displaystyle\vphantomu \over v\widehat\Gamma_gh(\omega) can be written as

6?, Phase Coherence

It is important, however, to emphasise that the similarity measure computed using the Parsimonious Cross Correlation of Waveforms C_l(t) uses raw amplitudes of ground movements (e.g., wave velocities).  For events that have similar signals, a Direct Cross Correltion Algorithm works well, even when signals get buried in large amounts of noise.  But when signals have different energy levels, strong signals will dominate the Parsimonious Cross Correlation Measure.  This means that with the Direct Running Cross Correlation, it becomes difficult to determine the similarity of two signals just by comparing the amplitude of their Cross Correlation, particularly when there are preprocessing biases that affect the amplitudes.

Apart from similarity measures based on raw amplitudes, various other metrics have been formulated throughout the years, each with its own strengths and weaknesses.  Some are based on normalisation techniques; whilst others utilise newer representations of signals that eliminate the weaknesses of the Fourier Transform.  Particularly relevant are Phase Coherence Similarity Measures (e.g, Hilbert Transforms) and Time-Frequency Spectral Transforms (e.g., Short Time Fourier Transforms, S Transforms, Wavelet Transforms).  The advantage of Wavelet Transforms is that besides offering better representations of Local Scale Physical Phenomena, they can be computed using very efficient algorithms.  As a consequence of the flexible nature of Wavelet Transforms, they are ideal for use in Ambient Based Microseismic Studies, particularly for Subsurface Tomography.

Schimmel (1999)Schimmel (1999): Phase Cross-Correlations: Design, comparisons, and applications.

The Phase Metric R_fg^\rm P(t') was motivated by the Phase Weighting Factors W^\rm P(t) introduced in the computation of the Phase Weighted Stack \overlineR^\rm P(t) to down-weight incoherent noise when computing Ordinary Direct Linear Stacks \overlineR^\rm L(t) (e.g., in the computation of Slant Stacks) in Schimmel & Paulssen (1997)Schimmel & Paulssen (1997): Noise reduction and detection of weak, coherent signals through phase weighted stacks..

6.7, Phase Coherency

igc06--PhaseWtStk show the Phase Weighted Stack described in Schimmel & Paulssen (1997)Schimmel & Paulssen (1997): Noise reduction and detection of weak, coherent signals through phase weighted stacks. applied to a Low Fold Coverage Seismic Profile, the different reflection arrivals identified by positive integers, both with the original synthetic (Left Vertical Panel) and with Background Noise (Right Vertical Panel).  The test results obtained by the Phase Weighted Stack are shown in the Lower Two Horizontal Panels.  The Central Horizontal Panel shows a stacking comparison along the Zero Slope Trajectory between: 1) the trace generated by the Phase Weighted Stacking; and 2) the trace obtained by using the Parsimonous Slant Stack.  The results for the Full Multi-Dip Version of the stacking algorithm can be seen at the Bottom Horizontal Panel.  Naturally, both the original synthetic and the noisy record are considered (simply compare the Left Vertical Panel with the Right Vertical Panel).  igc06--PhaseWtStk clearly shows that trace stacking using Local Phase Coherence is more effective in suppressing random noise and artifacts than other stacking methods, such as the Direct Linear Stack -  for both the original record and the noisy version.

6.8, Coherency on Quiet Signals

igc06--PhaseCoh--WvSignal shows three Coherence Functions between a Reference Signal and a Representative Seismic Trace.  The first function was determined using the Phase Coherence Algorithm; whilst the other two employed the Parsimonious Cross Correlation Method computed using the Normalised and Non-Normalised versions.  The Parsimonious Cross Correlation shows the strong sensitivity to large amplitude fluctuations in the Representative Trace.  But when the Phase Coherence is determined, only the correlation that exists between the phases is measured.  Consequently, the amplitude does not govern the heights of the peaks in the Phase Coherence Function.

6.9, Coherency of Agitated Waveforms

igc06--PhaseCoh--NoisyWvfrm picks two different wavetrains T_r and T_b from one seismic trace u_1(t).  The two wavetrains are then correlated with a second trace u_2(t).  The Phase Coherence using Wavetrain T_r picks up the correct phase shift.  Unlike the Normalised Parsimonious Cross Correlation - which fails - returning a phase shift of zero instead.

6.10, The Stockwell Transform

In 1996, the Stockwell Transform was introduced as a bridge between Short Time Fourier Transforms and Continuous Wavelet Transforms Stockwell, Mansinha & Lowe (1996)Stockwell, Mansinha & Lowe (1996): Localization of the complex spectrum: The S Transform..  Adopting shorter time windows for higher frequencies and larger ones for lower frequencies, the Stockwell Transform offers a better way to perform Fourier Based Frequency Analysis locally - analyzing dynamic frequency changes and a wider range of frequency resolutions than the Short Time Fourier Transform whilst maintaining Zero Phase at Zero Time.

In seismology the first signals to arrive from an earthquake are Body Waves - where P Waves and S Waves are followed by other Body Waves traveling along different paths and by higher amplitude Dispersive Surface Waves.  Because the spectral components clearly have a strong dependence on time, a Joint Time Frequency Representation (Tfr) is appropriate and unavoidable. (Joint in parameters t' and \omega.) For one dimensional signals, the Stockwell Transform is defined according to the following definition.

6.1, The Stockwell Transform

Definition

10a, The Gaussian and Hyperbolic Stockwell Transforms

10b, The Generalised Stockwell Transform

The Stockwell Transform is easily extended with additional parameters \\kappa_j\ by generalising the Window Function from W_t(t',\omega) to the arbitrary one W_t(t',\omega,\\kappa_j\).

, The Generalised Stockwell Transform

Using the Generalised Stockwell Transform, the Gaussian Stockwell Transform can be extended to include an additional parameter \gamma_G.

In a similar manner, the five parameter Hyperbolic Stockwell Transform is defined by Pinnegar & Mansinha (2003)Pinnegar & Mansinha (2003): The S-Transform with windows of arbitrary and varying shape..

6.11, Signal Recognition using Radon Transforms

The Stockwell Transform has recently been employed to correct the complex Doppler Frequency and Range Migration Problems associated with fast moving Hypersonic Weapons that present themselves even on short time scales Sun et al. (2019)Sun et al. (2019): Radon-S transform for hypersonic maneuvering target detection..  Hypersonic Weapons deliver conventional or nuclear payloads at ultra high velocities exceeding Mach 5 over long distance ranges, and are specifically designed to ameliorate their survivability against Ballistic Missile Defense Systems.  The novel target imaging technique is particularly useful against Highly Maneuvering Radar Targets in low Signal-To-Noise Ratio conditions.  The method uses the Stockwell Transform in conjuction with the traditional Radon Transform Method to update target lock on information obtained from Mid Course Guidance Signals in order to instruct the Missile Control System (e.g., Vane Control Systems, Thrust Vectoring Systems) to adjust the flight path.

The Radon Transform finds use in Multiple Reverberation Suppression and Signal Enhancement of seismic reflection data Shankar et al. (2009)Shankar et al. (2009): Signal enhancement and multiple suppression using Radon transform: An application to marine multichannel seismic data.  For Water Multiple Attenuation, the Parabolic Radon Transform has been widely used to integrate - in the mathematical sense of integration - Moveout Correction Common Midpoint Gathers (Cmp). Although the Parabolic Radon Transform is still quite unfamiliar to practitioners outside Exploration Geophysics, the Linear Path Version, popularly called the Slant Stack, has played a very important role in attenuating Ground Roll and in improving the Signal to Noise Ratio of seismic data Zhou & Greenhalgh (1994)Zhou & Greenhalgh (1994): Linear and parabolic t-p transforms revisited..

The Transforms have also dealt with the ill-posed and very hard to solve problem of Seismic Data Regularization, whereby seismic entire seismic traces are interpolated or extrapolation to regular grids - a requirement of many processing workflows used in Exploration Geophysics to infer the presence or position of economically useful geological deposits (e.g., fossil fuels and other hydrocarbons) Trad, Ulrych & Sacchi (2002)Trad, Ulrych & Sacchi (2002): Accurate interpolation with high-resolution time-variant Radon Transform.  The two Radon Transform Optimisers are customarily implemented in the Temporal Fourier Domain for faster computations, utilising the Toeplitz Structure coming from the Normal Equations in the Least Squares Inverse Estimator Yilmaz (2001)Yilmaz (2001): Seismic data analysis..

11a, Common Midpoint Stacking

In RadonTfm, it was argued that stacking Common Midpoint Gathers via the Radon Transform improves the Signal to Noise Ratio of seismic data by attenuating Incoherent Noise. (Common Midpoint Gathers are composed of traces that share the same point on the surface, halfway between source and receiver.)  In its simplest case, Common Midpoint Gathers are first Normal Moveout Corrected (Nmo Corrected), and then averaged horizontally.  Horizontal Stacking across the Gather forms a single Zero Offset Seismic Trace, reducing the amount of traces by a factor called the Fold of the Gather.

The main difficulties in stacking are: 1) All traces in a gather have equal validity; 2) All signals are coherent; and 3) All noises are random.  These assumptions, however, are valid only when: a) Cmp Sorting results in Midpoint Gathers having perfectly aligned seismic events from one trace to the next; b) Amplitudes are uniformly distributed across gather traces; and c) All noises are uncorrelated and stationary.  For all these reasons, the departure from the ideal Cmp Gather is well known but regularly ignored.  For large offsets, for instance, Hyperbolic Moveout (i.e., Normal Moveout of the actual traveltime with trace offset).  Moreover, the determination of Cmp Stacking Velocities is based on manually selecting a small number of locations (both Temporal [on Event Location] and Spatial [on Cmp Gather Location]) and picking Cmp Velocities that yields the best stack response at the selected events.  Using the Selected Cmp Gather Velocities as reference, the manually picked Cmp Velocities are extended for the whole survey area by interpolation (usually linear) between smoothed velocity curves.  The quality of the Cmp Stacks is therefore limited in good portion by human judgment that could be inadequate, particularly for large data sets. Another common problem that is encountered is trace stretching associated with Moveout Correction.  Thusly, traces that suffer from severe Moveout Stretch have their earlier amplitudes discared (i.e., zeroed out) Miller (1992)Miller (1992): Normal moveout stretch mute on shallow-reflection data. - introducing several problems such as: 1) Decreasing Fold, 2) Gaps in Three Dimensional Volumes; 3) Reduced Noise Suppression; 4) and Signal Destruction.

Consequently, there exists a need for subsurface imaging techniques that go beyond the rudimentary procedure of Parsimonious Stacking, in order to ameliorate quality and improve signal resolution.  This transition towards Non Parsimonious Seismic Processing is especially relevant to unconventional data types such as Multi Component Seismic Data, Extra Terrestrial Helioseismic and Asteroseismic Data.  And specifically for Coherence Data obtained from Microseismic Intersation Similarity Measures using Sophisticaled Microseismic Extraction Techniques based on: 1) Phase Coherence; 2) Stockwell Transforms; and 3) Wavelet Transforms.  These three special processing routines will surely benefit from unconventional stacking techniques.  In Exploration Geophysics, there is a wide class of Non Parsimonious Noise Suppression Techniques, the most interesting being Deconvolution, Migration Chun & Jacewitz (1981)Chun & Jacewitz (1981): Fundamentals of frequency-domain migration. (after the development on Kirchhoff Summation and Finite Difference Migration, Migration by Fourier Transforms was shown by Stolt (1978)Stolt (1978): Migration by Fourier Transform.), Surface Related Multiple Elimination and Adaptive Subtraction Verschuur (2013)Verschuur (2013): Seismic multiple removal techniques..

There are two main reasons why the development of Stacking Techniques have endured for so long.  Firstly, because in variety of cases stacking is the only process capable of separating signal from noise.  And most importantly, because many major seismic data processing steps rely some way or another on the concept of stacking (e.g, Velocity Analysis, Seisimic Migration).

6.12, The Stockwell Phase Coherence Stack

A Phase Coherence Stacking Approach based on the Stockwell Transform that characterises better the frequency components of Microseismic Signals was first proposed by Schimmel & Paulssen (1997)Schimmel & Paulssen (1997): Noise reduction and detection of weak, coherent signals through phase weighted stacks.  The suppression of incoherent microseismic noise from Daily Parsimonious Cross Correlation Records and from Daily Phase Coherence Records of Ambient Microseismicity is achieved by pointwise multiplication of a Linear Stack \overlineS(t',\omega) with a Time Dependent Weighting Filter W^G_t(t',\omega) that reduces the contribution of incoherent waveform features in the Stockwell Transform Domain.

The Time Averaged Stack \overlineR_S'(t) can then be obtained by an Inverse Stockwell Transform.

But, although one might get the idea that the Non Parsimonious Time Averaged Stack using the Stockwell Weighting Filter has been adequately realised, it is not so.  This is because Manipulated Stockwell Spectra (such as S\\overlineR_S'(t)) imposed by the Weighting Filter
W^G_t(t',\omega) suffer from severe artifacts after the Inverse
Stockwell Transformation is applied.
Schimmel & Gallart (2005)Schimmel & Gallart (2005): The inverse S Transform in
filters with time-frequency localization. have covered the inverse
mapping problem when using the routine Inverse Stockwell Transform by
Stockwell et al. (1996)Stockwell et al. (1996): Localization of the complex
spectrum: The S Transform., and propose a new approach that applies a
Frequency Dependent Scaling Factor k'(\omega) to the Stockwell
Spectrum.

Translating to the Angular Frequency \omega, the expression for S^-1\S(t,\omega)\ is rewritten as follows.

Examples with Synthetic Data and Real Data using Ocean Bottom Seismometers are discussed in Schimmel & Gallart (2007)Schimmel & Gallart (2007): Frequency-dependent phase coherence for noise suppression in seismic array data.  The Real Data Examples used: 1) A Combined Marine and Land Seismic Survey acquired by the 2003 Spanish Marconi Experiment Fernandez, Gallastegui & Pulgar (2009)Fernández, Gallastegui & Pulgar (2009): Structure and tectonic evolution of the eastern Cantabrian margin: Results from the Marconi multichannel seismic data. and 2) The May 26, 2003, M_W7 Earthquake near the East Coast of Honshu (Japan) of depth of 71.2km. (Marconi Experiment stands for North Iberian Continental Experiment (Spanish: Margen Continental Nor Iberico).  The Experiment consisted of: 1) A Wide Angle Hybrid Reflection Refraction Survey in the Bay of Biscay using Ocean Bottom Seismometers and Air Gun Sources [2.4km Linear Multi Source Shot Towed Array Streamer consisting of 96 Air Guns]) in the ; and 2) A Land Survey in the Cantabrian Mountain Range and Basque Pyrenees of Northern Spain.)

6.14, The Stockwell and Wavelet Connection

According to Ventosa et al. (2008)Ventosa et al. (2008): The S-Transform from a wavelet point of view.

14a, The Stockwell and Wavelet Transforms

6.2, The Wavelet Transform

where the family of waveforms \\psi_\alpha\beta(t)\ is generated by a Temporal Translational Parameter \alpha and a Temporal Scaling Parameter \beta of a Mother Wavelet \psi_M.

Definition

14b, The Formalisms of Fourier Transforms

Several formalisms exists for the Fourier Transform.

Fourier Formalism

Equipped with the Inner Product, the set of Fourier Basis \psi_k = \\cos \omega_k t, \sin \omega_k t \ are Orthogonal and Orthonormal.

A comparison can be made between the Fourier Transform with the Wavelet Transform by applying the notational formalisms of the Wavelet Transform to the Fourier Transform.  Herein, the Generalised form on the Fourier Transform is presented using the set of Wavelet Functions \ \psi_\omega(t) \.

14c, The Wavelet Based Stockwell Transform

Ventosa et al. (2008)Ventosa et al. (2008): The S-Transform from a wavelet point of view. present a case for improving the impractical computation costs involved in running the Phase Coherence Stack based on the Time Dependent Weighting Scheme constructed from the Gaussian Stockwell Transform of the Daily Interstation Phase Coherences.

The fundamental objection on using the Wavelet Transform for extracting Microseismic Wavetrains excited by Atmosphere-Ocean-SeaFloor Coupling is its diminishing amplitudes towards higher frequencies igc06--GsStkwl--MrltWvlt.  But, the authors, in their creative phase, rejuvenated the interest in Wavelet Transforms for Microseismology by explaining how the Morlet Wavelet Spectrum and the Stockwell Spectrum are related.  In particular, they derive the Two Weighting Matrices C_ws(t',\omega) and C_sw(t',\omega) that allow a deterministic mapping between the two Time Frequency Transforms - solving the Amplitude Variablity Problem when using the Wavelet Transform igc06--GsStkwl--MrltWvlt.

igc06--GsStkwl--MrltWvlt--NsyWvfrm shows an additional synthetic example - using a Chirplet contaminated with a Gaussian Probability Density where the contaminating values are statistically uncorrelated and identically distributed. (In seismology, Guassian Noise is used as additive White Noise that characteristically has a Constant Power Spectral Density).  The agreeement between the two amplitude spectra - One being the Gaussian Stockwell Spectrum; the other being the Modified Morlet Wavelet Transform mapped on frequency rather than scale - is unassailable.  As deduced by Ventosa et al. (2008)Ventosa et al. (2008): The S-Transform from a wavelet point of view., C_ws(t',\omega) produces the necessary Amplitude Correction to the Morlet Wavlet Transform so as to be useful for extracting embedded microseismic signals within ground motion seismograms.

6.15, Sicilian Microseismic Wavefield Extraction by
Morlet Wavelet Phase Coherence Stacking

igc06--Siculo--WvltCoh--Stk are the results obtained from computing the three Coherence Stacks discussed in this transfer report. The illustration starts with a progressive averaging scheme that estimates the Daily Parsimonious Cross Correlation Function computed on continuous vertical component seismograms selected from seismic stations in Sicily.

The Parsimonious Normalised Correlation Function is computed in the frequency domain by running a Fast Fourier Transform Algorithm, then reverting back to the time domain.  The Ambient Microseismic Pre-Processing Scheme between station seismograms was performed on the largest interstation distances as discussed in Chapter 03.  One Bit Normalisation on seismograms is used to reduce the contribution of the most energetic motions excited by earthquakes.

A comparison is made between the final results obtained by stacking the computed cross correlations functions for a set of continuous recording temporal periods - 2months, 4months, 6months, 8months, and 1year.  See igc06--Siculo--WvltCoh--Stk.  Comparisons between the five stacked traces for R_l(t) are used as quality control metrics to determine the averaging that is required to effectively enhance the Signal-to-Noise Ratio from one stacked trace to the next.  Through these comparisons, it is concluded that the meteorological phenomena observed over the Mediterranean Sea are mild.  Only on very rare occasions, do storms reach the strength of a Category 1 Hurricanes (119km/h to 153km/h) as defined by the Saffir Simpson Hurricane Wind Scale Schott et al. (2019).  The Asymmetry is consistent with the prevalent wind and wave conditions in the Mediterranean Sea Cavaleri (2005).  See igc06--MdSeaWvHs.

igc06--Siculo--WvltCoh--Stk also shows the results obtained for the Phase Coherence Stacks obtained by running the Gaussian Stockwell Transform and the Time-Scale Morlet Wavelet Transform, for One Year of Continuous Vertical Ground Motion Seismograms.  The enhanced Signal-to-Noise Level and the reduced noise components provides confirmation of how the measurements of the Ambient Wavefield Phase present a powerful tool for extracting weak surface wave microseismic wavetrains at high resolution.  The figure also shows the similarity between the two Time Frequency Representations and how the Inverse Stockwell Transform does indeed produce a Coherence-Like Trace.

Inspecting vertical component seismograms, stacked sicilian correlograms of ambient microseismicity are able to identify dominant noise sources near the North Western Mediterranean Sea (NW-MS) Deep Water Formation as a moderate source of sea storms Picotti et al. (2014).  A comparison of the Group Velocity obtained from the Parsimonious Cross Correlation Model and the Time-Scale Morlet Wavelet Transform is shown in igc06--SpectrDispr--CohStk.  The coherence derived Droup Velocity Dispersion Curve Fitting Results identified by the solid black curves give a good representation of the frequency dispersion characteristic of Earth Surface Rayleigh Wave Propagation. Furthermore, an interpretation of the Velocity-Frequency Density Pattern identifies a region of significant contribution from Fundamental Rayleigh Waves consistent with Secondary Microseisms located at the familiar 140mHz region.  The prevailing region of significant microseismic contributions driven by Atmosphere-Ocean-SeaFloor Coupling is the result of storms in the Western Mediterranean Sea.  The second conclusion on microseismicity in the Central Mediterranean is that the contributions from Omnidirectional Ambient Noise Measurements acquired in Sicily are for most part not associated with surf generated noise characteristic of Primary Microseisms.

Nevertheless, it is vital to grasp the reality that satisfactory extracting Group and Phase Velocities from weak Microseismic Signals to understand the subsurface is extremely challenging.  It is therfore clear that an understanding of the weak and subtle ground motion observations in sync with discussions on the hottest topics (outlined in HotShit--Intgr--SmPract) discussed and fully debated.

6.16, The Hot Topics for an
Integrative Approach to Seismological Practice

The improvements - both seismological and computational -to seismological practice, using as example the analysis of continuous ground reverberations for the extraction and use of weak microseismic signals excited by Atmosphere-Ocean-SeaFloor Coupling, illustrates the value of an integrative approach spanings the continuum of understanding from fundamental considerations to the practical setup that enhances understanding and free collaboration.  To achieve this goal, five hot topics require specialised elucidation.

Thusly, I present my particular view that is concentrated on using large continuous seismological datasets for analysing and applying the weak background microseismic wavetrains that are customarily hidden from view because of the noisy character of seismograms recorded by the global seismological networks in operation.  I currently offer my specific view on how Seismological Practice is today, but my discussion is structured in ways that will remain contemporary with the evolving and emerging fields within seismology and beyond.  Still, it remains difficult to advance a deep and complete analysis of the relevant discussions that are so rich in the variety of topics that would have to be covered.

Because of the timeframes imposed for the completion of this report, the five sections shall be completed afterwards.  Nevertheless, this transfer report calls upon the expertise of an extensive list of professionals across multiple colleges to an integrating apporach to become standard practice.

16A, A focus on Tools that Respect the Freedom of Users

16B, The Application of Mathematical Intuition in contrast to Purely Physical Intuitions to answer Current and Future Seismological Questions

16C, A Drive towards Automation of Computational Pipelines

16D, The Utilisation and Free Sharing of Comprehensive Sets of Rapid High End Computing Techniques

16E, An Integrative Organisational Framework to Seismological Practice that drives Free Collaboration of Value

6B IDEOGRAMS CITED IN CHAPTER 6
*******************************

6C REFERENCES CITED IN CHAPTER 6
********************************

6D CREATIVE COMMONS (BY-NC-ND)
******************************

4cm Creative Commons
1cm Attribution-NonCommercial-NoDerivatives
3.5cm 4.0 International License