Visualization and Imaging

# Trace plotting methods
Trace.plot(self, type: str = 'normal', starttime=None, endtime=None, 
          fig=None, ax=None, title: str = None, show: bool = True, 
          draw: bool = None, block: bool = None, **kwargs):
    """
    Plot single trace waveform.
    
    Args:
        type: Plot type ('normal', 'dayplot', 'relative')
        starttime: Start time for plot window (UTCDateTime or None)
        endtime: End time for plot window (UTCDateTime or None)
        fig: Matplotlib figure object (created if None)
        ax: Matplotlib axes object (created if None) 
        title: Plot title (auto-generated if None)
        show: Display plot immediately
        draw: Force matplotlib draw
        block: Block execution until window closed
        **kwargs: Additional matplotlib plotting options
        
    Matplotlib kwargs:
        color, linewidth, linestyle, alpha, label, etc.
    """

# Stream plotting methods  
Stream.plot(self, type: str = 'normal', starttime=None, endtime=None,
           fig=None, size: tuple = (800, 250), title: str = None,
           color: str = 'black', number_of_ticks: int = 4, 
           tick_rotation: int = 0, tick_format: str = None,
           show: bool = True, draw: bool = None, block: bool = None,
           equal_scale: bool = True, **kwargs):
    """
    Plot multiple traces in stream.
    
    Args:
        type: Plot type ('normal', 'dayplot', 'section')
        starttime: Plot window start time
        endtime: Plot window end time
        fig: Matplotlib figure object
        size: Figure size in pixels (width, height)
        title: Plot title
        color: Trace color (single color or list)
        number_of_ticks: Number of time axis ticks
        tick_rotation: Rotation angle for tick labels
        tick_format: Time tick format string
        show: Display plot
        draw: Force draw
        block: Block until closed
        equal_scale: Use same amplitude scale for all traces
        **kwargs: Additional plotting options
        
    Plot Types:
        'normal': Standard time series plot
        'dayplot': 24-hour record section plot  
        'section': Record section with distance/time axes
    """

# Trace spectrogram
Trace.spectrogram(self, log: bool = False, outfile: str = None, 
                 format: str = None, axes=None, show: bool = True,
                 title: str = None, **kwargs):
    """
    Plot spectrogram of trace data.
    
    Args:
        log: Use logarithmic frequency scale
        outfile: Save plot to file (format from extension)
        format: Image format for saving
        axes: Matplotlib axes object
        show: Display plot
        title: Plot title
        **kwargs: Spectrogram parameters
        
    Spectrogram Parameters:
        samp_rate: Sampling rate override
        per_lap: Window overlap percentage (0-1)
        wlen: Window length in seconds
        mult: Frequency scaling multiplier
        cmap: Colormap name
        vmin, vmax: Color scale limits
        dbscale: Use decibel scale
    """

# Stream spectrogram
Stream.spectrogram(self, **kwargs):
    """
    Plot spectrograms for all traces in stream.
    
    Args:
        **kwargs: Same as Trace.spectrogram()
        
    Creates subplot for each trace in stream.
    """

# Import from obspy.imaging.beachball
def beachball(fm, linewidth: int = 2, facecolor: str = 'b', bgcolor: str = 'w',
             edgecolor: str = 'k', alpha: float = 1.0, xy: tuple = (0, 0),
             width: int = 200, size: int = 100, nofill: bool = False,
             zorder: int = 100, outfile: str = None, format: str = None,
             fig=None, **kwargs):
    """
    Draw focal mechanism beachball.
    
    Args:
        fm: Focal mechanism (strike, dip, rake) or moment tensor
        linewidth: Nodal plane line width
        facecolor: Compressional quadrant color
        bgcolor: Background color  
        edgecolor: Outline color
        alpha: Transparency (0-1)
        xy: Center position for plot
        width: Beachball width in points
        size: Beachball size in points
        nofill: Draw outline only
        zorder: Matplotlib drawing order
        outfile: Save to file
        format: Image format
        fig: Matplotlib figure
        **kwargs: Additional options
        
    Focal Mechanism Formats:
        [strike, dip, rake]: Single nodal plane (degrees)
        [strike1, dip1, rake1, strike2, dip2, rake2]: Both planes
        MomentTensor object: Full moment tensor
        [mrr, mtt, mpp, mrt, mrp, mtp]: Moment tensor components
    """

def beach(fm, **kwargs):
    """
    Create beachball plot (alias for beachball function).
    
    Args:
        fm: Focal mechanism data
        **kwargs: Same as beachball()
        
    Returns:
        Matplotlib collection object
    """

# Create day plot from continuous data
st = read('continuous_24h.mseed')
st.plot(type='dayplot', 
        starttime=UTCDateTime("2023-01-01T00:00:00"),
        title="24-Hour Seismic Record",
        color=['red', 'blue', 'green'],  # Different colors per channel
        number_of_ticks=6,  # 4-hour intervals
        tick_format='%H:%M',
        size=(1200, 800))

# Record section plotting parameters for Stream.plot()
SECTION_PLOT_OPTIONS = {
    'type': 'section',
    'recordlength': 300,      # Trace length in seconds  
    'recordstart': 0,         # Start offset in seconds
    'norm_method': 'trace',   # Normalization ('trace', 'stream', 'global')
    'distance_degree': True,  # Use degree distances
    'ev_coord': tuple,        # Event coordinates (lat, lon)
    'alpha': 0.8,            # Transparency
    'linewidth': 1.0         # Line width
}

# PPSD plotting methods (see signal-processing.md for PPSD class)
PPSD.plot(self, filename: str = None, show_coverage: bool = True,
         show_histogram: bool = True, show_percentiles: bool = False, 
         percentiles: list = [0, 25, 50, 75, 100], show_noise_models: bool = True,
         grid: bool = True, period_lim: tuple = None, show_mean: bool = False,
         **kwargs):
    """
    Plot PPSD analysis results.
    
    Args:
        filename: Save plot to file
        show_coverage: Show temporal coverage
        show_histogram: Show probability density histogram  
        show_percentiles: Show percentile curves
        percentiles: List of percentiles to display
        show_noise_models: Show Peterson noise models
        grid: Show grid lines
        period_lim: Period axis limits (min, max)
        show_mean: Show mean power levels
        **kwargs: Additional plotting options
    """

PPSD.plot_temporal(self, starttime=None, endtime=None, filename: str = None,
                  **kwargs):
    """
    Plot temporal evolution of power spectral densities.
    
    Args:
        starttime: Time window start
        endtime: Time window end  
        filename: Save to file
        **kwargs: Plotting options
    """

PPSD.plot_spectrogram(self, filename: str = None, clim: tuple = None, **kwargs):
    """
    Plot spectrogram view of PPSD data.
    
    Args:
        filename: Output filename
        clim: Color scale limits
        **kwargs: Additional options
    """

from obspy import read, UTCDateTime
import matplotlib.pyplot as plt

# Read seismic data
st = read('seismic_data.mseed')

# Simple trace plot
trace = st[0]
trace.plot()

# Customized trace plot
trace.plot(color='red', linewidth=1.5, 
          title=f"Seismogram: {trace.id}",
          starttime=trace.stats.starttime,
          endtime=trace.stats.starttime + 300)  # First 5 minutes

# Stream plot with multiple traces
st.plot(type='normal', 
        size=(1000, 600),
        color=['black', 'red', 'blue'], 
        equal_scale=True,
        title="Three-Component Seismogram")

from obspy import read
import matplotlib.pyplot as plt

# Read earthquake data
st = read('earthquake_record.mseed')
trace = st.select(channel="*Z")[0]  # Vertical component

# Basic spectrogram
trace.spectrogram(show=True)

# Advanced spectrogram with custom parameters
fig, ax = plt.subplots(figsize=(12, 8))
trace.spectrogram(
    log=True,           # Logarithmic frequency scale
    axes=ax,
    show=False,
    # Spectrogram parameters
    wlen=10.0,          # 10-second windows
    per_lap=0.9,        # 90% overlap
    mult=2.0,           # Frequency multiplier
    cmap='plasma',      # Color scheme
    vmin=-160,          # Color scale minimum
    vmax=-100,          # Color scale maximum  
    dbscale=True        # Decibel scale
)
ax.set_title(f"Spectrogram: {trace.id}")
ax.set_ylabel("Frequency (Hz)")
ax.set_xlabel("Time")
plt.show()

# Multi-channel spectrograms
st_3c = st.select(station="ANMO")  # All components from one station
st_3c.spectrogram()  # Creates subplot for each component

from obspy.imaging.beachball import beachball, beach
import matplotlib.pyplot as plt

# Focal mechanism parameters (strike, dip, rake)
focal_mechanisms = [
    [95, 85, 180],      # Strike-slip
    [0, 90, -90],       # Normal fault  
    [45, 30, 90],       # Reverse fault
    [315, 45, 45],      # Oblique
]

mechanism_types = ["Strike-slip", "Normal", "Reverse", "Oblique"]

fig, axes = plt.subplots(2, 2, figsize=(10, 10))
axes = axes.flatten()

for i, (fm, mtype) in enumerate(zip(focal_mechanisms, mechanism_types)):
    # Create beachball on specific axes
    ax = axes[i]
    beachball(fm, ax=ax, size=100, 
             facecolor='red',     # Compressional quadrants
             bgcolor='white',     # Background
             edgecolor='black',   # Outline
             linewidth=2)
    
    ax.set_title(f"{mtype}\nStrike={fm[0]}°, Dip={fm[1]}°, Rake={fm[2]}°")
    ax.set_xlim(-1.2, 1.2)
    ax.set_ylim(-1.2, 1.2)
    ax.set_aspect('equal')

plt.tight_layout()
plt.show()

# Beachballs with moment tensor
from obspy.imaging.beachball import MomentTensor

# Harvard CMT solution (Mrr, Mtt, Mpp, Mrt, Mrp, Mtp)
mt_components = [1.0, -1.0, 0.0, 0.0, 0.0, 0.0]  # Simple strike-slip

beachball(mt_components, 
         outfile='focal_mechanism.png',
         format='PNG',
         size=200)

from obspy import read, UTCDateTime

# Read 24-hour continuous data
st = read('continuous_data_24h.mseed') 

# Create day plot
start_time = UTCDateTime("2023-01-15T00:00:00")
st.plot(type='dayplot',
        starttime=start_time,
        title="24-Hour Seismic Record - Station ANMO",
        size=(1400, 1000),
        color='black',
        number_of_ticks=12,  # 2-hour intervals
        tick_format='%H:%M',
        show_y_UTC_label=True)

# Day plot with spectrogram overlay
st_single = st.select(channel="BHZ")  # Single channel
st_single.plot(type='dayplot',
              starttime=start_time, 
              title="Vertical Component with Spectrogram",
              show_spectrogram=True,
              spectrogram_cmap='plasma')

from obspy import read, UTCDateTime
from obspy.geodetics import gps2dist_azimuth

# Read multi-station earthquake data
st = read('regional_earthquake.mseed')

# Event location
event_lat, event_lon = 34.0, -118.0
event_time = UTCDateTime("2023-01-15T10:30:00")

# Calculate distances and add to trace stats
for trace in st:
    # Get station coordinates (from StationXML or manual)
    sta_lat = trace.stats.coordinates.latitude  
    sta_lon = trace.stats.coordinates.longitude
    
    # Calculate distance
    dist_m, az, baz = gps2dist_azimuth(event_lat, event_lon, sta_lat, sta_lon)
    trace.stats.distance = dist_m / 1000.0  # km

# Sort by distance
st.sort(['distance'])

# Create record section
st.plot(type='section',
        recordlength=180,        # 3-minute traces
        recordstart=-30,         # Start 30s before origin
        norm_method='trace',     # Normalize each trace
        distance_degree=False,   # Use km distances
        ev_coord=(event_lat, event_lon),
        alpha=0.8,
        linewidth=1.0,
        title="Regional Earthquake Record Section")

from obspy import read
from obspy.signal import PPSD
from obspy.clients.fdsn import Client

# Get instrument response
client = Client("IRIS")
inventory = client.get_stations(network="IU", station="ANMO", 
                               location="00", channel="BHZ",
                               level="response")

# Read continuous data
st = read('noise_data.mseed')

# Create PPSD analysis
ppsd = PPSD(st[0].stats, metadata=inventory)
ppsd.add(st)

# Standard PPSD plot
ppsd.plot(show_coverage=True,
         show_histogram=True, 
         show_percentiles=True,
         percentiles=[10, 50, 90],
         show_noise_models=True,
         filename='ppsd_analysis.png')

# Temporal evolution plot
ppsd.plot_temporal(filename='ppsd_temporal.png')

# Spectrogram view  
ppsd.plot_spectrogram(filename='ppsd_spectrogram.png')

# Custom PPSD plot with multiple percentiles
import matplotlib.pyplot as plt

fig, ax = plt.subplots(figsize=(12, 8))
ppsd.plot(axes=ax, show=False,
         show_percentiles=True,
         percentiles=[5, 25, 50, 75, 95],
         grid=True,
         period_lim=(0.01, 100))

ax.set_title("Long-term Noise Analysis - Station ANMO")
ax.set_xlabel("Period (s)")
ax.set_ylabel("Power Spectral Density (dB)")
plt.show()

from obspy import read, UTCDateTime
import matplotlib.pyplot as plt
import matplotlib.dates as mdates

# Read high-quality data
st = read('publication_data.mseed')
trace = st[0]

# Create publication figure
fig, (ax1, ax2) = plt.subplots(2, 1, figsize=(12, 10))

# Waveform plot
times = trace.times("matplotlib")
ax1.plot(times, trace.data, 'k-', linewidth=0.8)
ax1.set_ylabel('Amplitude (counts)', fontsize=12)
ax1.set_title(f'Seismogram: {trace.id}', fontsize=14, fontweight='bold')
ax1.grid(True, alpha=0.3)

# Format time axis
ax1.xaxis.set_major_formatter(mdates.DateFormatter('%H:%M:%S'))
ax1.xaxis.set_major_locator(mdates.MinuteLocator(interval=2))
plt.setp(ax1.xaxis.get_majorticklabels(), rotation=45)

# Spectrogram
trace.spectrogram(axes=ax2, show=False,
                 wlen=30.0, per_lap=0.9,
                 log=True, dbscale=True,
                 cmap='viridis')
ax2.set_ylabel('Frequency (Hz)', fontsize=12)
ax2.set_xlabel('Time (UTC)', fontsize=12)

# Adjust layout and save
plt.tight_layout()
plt.savefig('publication_figure.pdf', dpi=300, bbox_inches='tight')
plt.savefig('publication_figure.png', dpi=300, bbox_inches='tight')
plt.show()

# Plot configuration structure
PlotConfig = {
    'type': str,              # Plot type ('normal', 'dayplot', 'section')
    'size': tuple[int, int],  # Figure size (width, height) pixels
    'color': str | list[str], # Colors for traces
    'linewidth': float,       # Line width
    'title': str,             # Plot title
    'starttime': UTCDateTime,  # Time window start
    'endtime': UTCDateTime,    # Time window end
    'equal_scale': bool,      # Equal amplitude scaling
    'number_of_ticks': int,   # Time axis tick count
    'tick_format': str        # Time tick format string
}

# Spectrogram parameters
SpectrogramConfig = {
    'wlen': float,            # Window length in seconds
    'per_lap': float,         # Overlap fraction (0-1)
    'mult': float,            # Frequency multiplier
    'cmap': str,              # Matplotlib colormap
    'vmin': float,            # Color scale minimum
    'vmax': float,            # Color scale maximum
    'log': bool,              # Logarithmic frequency scale
    'dbscale': bool           # Decibel amplitude scale
}

# Beachball styling
BeachballStyle = {
    'size': int,              # Beachball size in points
    'linewidth': float,       # Nodal plane line width
    'facecolor': str,         # Compressional quadrant color
    'bgcolor': str,           # Background color
    'edgecolor': str,         # Outline color
    'alpha': float,           # Transparency (0-1)
    'nofill': bool           # Outline only
}