Miserlou · April 22, 2026 17:48
diff --git a/batchplot.sh b/batchplot.sh
 #!/bin/bash

 TARGET_DIR="${1:-.}"
 SIZE_LIMIT=$((10 * 1024 * 1024))  # 10MB in bytes

 find "$TARGET_DIR" -type f | while read -r file; do
    # Get file size in bytes (works on OSX)
    filesize=$(stat -f%z "$file" 2>/dev/null)
    if [ -n "$filesize" ] && [ "$filesize" -gt "$SIZE_LIMIT" ]; then
        echo "Processing: $file ($filesize bytes)"
        python plot.py "$file"
    fi
 done
diff --git a/infraplot.py b/infraplot.py
 #!/usr/bin/env python3
 """
 infrasound_separate.py - Output two PNGs with timestamp in filename.
 Usage: python infrasound_separate.py input_file [output_directory]

 Slop..

 """

 import sys
 import os
 import numpy as np
 import matplotlib
 matplotlib.use('Agg')
 import matplotlib.pyplot as plt
 import matplotlib.dates as mdates
 from datetime import datetime, timedelta
 from obspy import read
 from scipy.signal import spectrogram

 def sanitize_timestamp(dt):
    """Convert datetime to safe filename string: YYYYMMDD_HHMMSS"""
    return dt.strftime('%Y%m%d_%H%M%S')

 def main():
    if len(sys.argv) < 2:
        print("Usage: python infrasound_separate.py <input_file> [output_directory]")
        sys.exit(1)

    infile = sys.argv[1]
    output_dir = sys.argv[2] if len(sys.argv) > 2 else os.path.dirname(infile)

    if output_dir and not os.path.exists(output_dir):
        os.makedirs(output_dir)

    # Read data
    st = read(infile)
    tr = st[0]
    data = tr.data.astype(np.float64)
    fs = tr.stats.sampling_rate
    starttime = tr.stats.starttime.datetime

    # Create filename with timestamp
    station_id = os.path.basename(infile).split('.')[0:3]  # e.g., ['AM', 'R80B3', '00']
    station_str = '.'.join(station_id) if station_id else 'unknown'
    timestamp_str = sanitize_timestamp(starttime)

    base_name = f"{station_str}_{timestamp_str}"

    print(f"File: {infile}")
    print(f"Station: {station_str}")
    print(f"Start time: {starttime}")
    print(f"Samples: {len(data)}, Sample rate: {fs} Hz, Duration: {len(data)/fs/3600:.2f} hours")

    # Time axis
    time_seconds = np.arange(len(data)) / fs
    time_datetime = [starttime + timedelta(seconds=t) for t in time_seconds]
    time_num = mdates.date2num(time_datetime)

    # ---------- 1. Waveform only ----------
    fig1, ax1 = plt.subplots(figsize=(14, 4))
    ax1.plot(time_num, data, 'b-', linewidth=0.5)
    ax1.set_ylabel('Amplitude (counts)')
    ax1.set_title(f'Infrasound Waveform - {station_str}\nStart: {starttime}')
    ax1.grid(True, alpha=0.3)

    # Format x-axis
    ax1.xaxis.set_major_formatter(mdates.DateFormatter('%H:%M:%S'))
    ax1.xaxis.set_major_locator(mdates.AutoDateLocator(maxticks=8))
    ax1.set_xlabel('Time (UTC)')
    plt.setp(ax1.xaxis.get_majorticklabels(), rotation=30, ha='right')

    plt.tight_layout()
    waveform_out = os.path.join(output_dir, f"{base_name}_waveform.png")
    plt.savefig(waveform_out, dpi=150, bbox_inches='tight')
    plt.close()
    print(f"Saved: {waveform_out}")

    # ---------- 2. Spectrogram only (0–10 Hz) ----------
    nperseg = min(1024, len(data)//8)
    f, t, Sxx = spectrogram(data, fs, nperseg=nperseg, noverlap=nperseg//2,
                            mode='psd', scaling='density')
    Sxx_dB = 10 * np.log10(Sxx + 1e-12)

    # Convert t to datetime
    t_datetime = [starttime + timedelta(seconds=ti) for ti in t]
    t_num = mdates.date2num(t_datetime)

    # Limit frequency to 0–10 Hz
    freq_mask = f <= 10
    f_plot = f[freq_mask]
    Sxx_plot = Sxx_dB[freq_mask, :]

    fig2, ax2 = plt.subplots(figsize=(14, 5))
    im = ax2.pcolormesh(t_num, f_plot, Sxx_plot, shading='gouraud', cmap='viridis')
    ax2.set_ylabel('Frequency (Hz)')
    ax2.set_title(f'Spectrogram (0–10 Hz) - {station_str}\nStart: {starttime}')
    ax2.set_ylim(0, 10)
    plt.colorbar(im, ax=ax2, label='Power (dB/Hz)')

    # Format x-axis
    ax2.xaxis.set_major_formatter(mdates.DateFormatter('%H:%M:%S'))
    ax2.xaxis.set_major_locator(mdates.AutoDateLocator(maxticks=16))
    ax2.set_xlabel('Time (UTC)')
    plt.setp(ax2.xaxis.get_majorticklabels(), rotation=30, ha='right')

    plt.tight_layout()
    spec_out = os.path.join(output_dir, f"{base_name}_spectrogram.png")
    plt.savefig(spec_out, dpi=150, bbox_inches='tight')
    plt.close()
    print(f"Saved: {spec_out}")

 if __name__ == '__main__':
    main()
diff --git a/requirements.txt b/requirements.txt
 certifi==2026.2.25
 charset-normalizer==3.4.7
 contourpy==1.3.3
 cycler==0.12.1
 decorator==5.2.1
 fonttools==4.62.1
 greenlet==3.3.2
 h5py==3.16.0
 idna==3.11
 kiwisolver==1.5.0
 lxml==6.0.2
 matplotlib==3.10.8
 numpy==2.4.4
 obspy==1.5.0
 packaging==26.0
 pillow==12.2.0
 pyparsing==3.3.2
 python-dateutil==2.9.0.post0
 requests==2.33.1
 scipy==1.17.1
 six==1.17.0
 SQLAlchemy==2.0.49
 typing_extensions==4.15.0
 urllib3==2.6.3
	#!/bin/bash

	TARGET_DIR="${1:-.}"
	SIZE_LIMIT=$((10 * 1024 * 1024)) # 10MB in bytes

	find "$TARGET_DIR" -type f \| while read -r file; do
	# Get file size in bytes (works on OSX)
	filesize=$(stat -f%z "$file" 2>/dev/null)
	if [ -n "$filesize" ] && [ "$filesize" -gt "$SIZE_LIMIT" ]; then
	echo "Processing: $file ($filesize bytes)"
	python plot.py "$file"
	fi
	done
	#!/usr/bin/env python3
	"""
	infrasound_separate.py - Output two PNGs with timestamp in filename.
	Usage: python infrasound_separate.py input_file [output_directory]

	Slop..

	"""

	import sys
	import os
	import numpy as np
	import matplotlib
	matplotlib.use('Agg')
	import matplotlib.pyplot as plt
	import matplotlib.dates as mdates
	from datetime import datetime, timedelta
	from obspy import read
	from scipy.signal import spectrogram

	def sanitize_timestamp(dt):
	"""Convert datetime to safe filename string: YYYYMMDD_HHMMSS"""
	return dt.strftime('%Y%m%d_%H%M%S')

	def main():
	if len(sys.argv) < 2:
	print("Usage: python infrasound_separate.py <input_file> [output_directory]")
	sys.exit(1)

	infile = sys.argv[1]
	output_dir = sys.argv[2] if len(sys.argv) > 2 else os.path.dirname(infile)

	if output_dir and not os.path.exists(output_dir):
	os.makedirs(output_dir)

	# Read data
	st = read(infile)
	tr = st[0]
	data = tr.data.astype(np.float64)
	fs = tr.stats.sampling_rate
	starttime = tr.stats.starttime.datetime

	# Create filename with timestamp
	station_id = os.path.basename(infile).split('.')[0:3] # e.g., ['AM', 'R80B3', '00']
	station_str = '.'.join(station_id) if station_id else 'unknown'
	timestamp_str = sanitize_timestamp(starttime)

	base_name = f"{station_str}_{timestamp_str}"

	print(f"File: {infile}")
	print(f"Station: {station_str}")
	print(f"Start time: {starttime}")
	print(f"Samples: {len(data)}, Sample rate: {fs} Hz, Duration: {len(data)/fs/3600:.2f} hours")

	# Time axis
	time_seconds = np.arange(len(data)) / fs
	time_datetime = [starttime + timedelta(seconds=t) for t in time_seconds]
	time_num = mdates.date2num(time_datetime)

	# ---------- 1. Waveform only ----------
	fig1, ax1 = plt.subplots(figsize=(14, 4))
	ax1.plot(time_num, data, 'b-', linewidth=0.5)
	ax1.set_ylabel('Amplitude (counts)')
	ax1.set_title(f'Infrasound Waveform - {station_str}\nStart: {starttime}')
	ax1.grid(True, alpha=0.3)

	# Format x-axis
	ax1.xaxis.set_major_formatter(mdates.DateFormatter('%H:%M:%S'))
	ax1.xaxis.set_major_locator(mdates.AutoDateLocator(maxticks=8))
	ax1.set_xlabel('Time (UTC)')
	plt.setp(ax1.xaxis.get_majorticklabels(), rotation=30, ha='right')

	plt.tight_layout()
	waveform_out = os.path.join(output_dir, f"{base_name}_waveform.png")
	plt.savefig(waveform_out, dpi=150, bbox_inches='tight')
	plt.close()
	print(f"Saved: {waveform_out}")

	# ---------- 2. Spectrogram only (0–10 Hz) ----------
	nperseg = min(1024, len(data)//8)
	f, t, Sxx = spectrogram(data, fs, nperseg=nperseg, noverlap=nperseg//2,
	mode='psd', scaling='density')
	Sxx_dB = 10 * np.log10(Sxx + 1e-12)

	# Convert t to datetime
	t_datetime = [starttime + timedelta(seconds=ti) for ti in t]
	t_num = mdates.date2num(t_datetime)

	# Limit frequency to 0–10 Hz
	freq_mask = f <= 10
	f_plot = f[freq_mask]
	Sxx_plot = Sxx_dB[freq_mask, :]

	fig2, ax2 = plt.subplots(figsize=(14, 5))
	im = ax2.pcolormesh(t_num, f_plot, Sxx_plot, shading='gouraud', cmap='viridis')
	ax2.set_ylabel('Frequency (Hz)')
	ax2.set_title(f'Spectrogram (0–10 Hz) - {station_str}\nStart: {starttime}')
	ax2.set_ylim(0, 10)
	plt.colorbar(im, ax=ax2, label='Power (dB/Hz)')

	# Format x-axis
	ax2.xaxis.set_major_formatter(mdates.DateFormatter('%H:%M:%S'))
	ax2.xaxis.set_major_locator(mdates.AutoDateLocator(maxticks=16))
	ax2.set_xlabel('Time (UTC)')
	plt.setp(ax2.xaxis.get_majorticklabels(), rotation=30, ha='right')

	plt.tight_layout()
	spec_out = os.path.join(output_dir, f"{base_name}_spectrogram.png")
	plt.savefig(spec_out, dpi=150, bbox_inches='tight')
	plt.close()
	print(f"Saved: {spec_out}")

	if __name__ == '__main__':
	main()
	certifi==2026.2.25
	charset-normalizer==3.4.7
	contourpy==1.3.3
	cycler==0.12.1
	decorator==5.2.1
	fonttools==4.62.1
	greenlet==3.3.2
	h5py==3.16.0
	idna==3.11
	kiwisolver==1.5.0
	lxml==6.0.2
	matplotlib==3.10.8
	numpy==2.4.4
	obspy==1.5.0
	packaging==26.0
	pillow==12.2.0
	pyparsing==3.3.2
	python-dateutil==2.9.0.post0
	requests==2.33.1
	scipy==1.17.1
	six==1.17.0
	SQLAlchemy==2.0.49
	typing_extensions==4.15.0
	urllib3==2.6.3