Course 5 - Visualize earthquake peaks

clear all
close all
clc

% 1. Power in dB from [p, f, t] = [power, frequency, time]
load after_filtering
[p, f, t] = pspectrum(wanc_bandpass,50,"spectrogram")
psum = sum(p);
plot(t, psum)
pwr = db(psum, "power")
plot(t, pwr)

% 2. Local extrema live task - Find local maxima here
maxIndices = islocalmax(pwr,"MinProminence",10,"SamplePoints",t);
figure
plot(t,pwr,"Color",[77 190 238]/255,"DisplayName","Input data")
hold on
plot(t(maxIndices),pwr(maxIndices),"^","Color",[217 83 25]/255,...
    "MarkerFaceColor",[217 83 25]/255,"DisplayName","Local maxima")
title("Number of extrema: " + nnz(maxIndices))
hold off
legend
xlabel("t")

% 3. Visualize earthquakes
quakeSec = t(maxIndices)
tvec = (0:numel(wanc_bandpass)-1)/50;
plot(tvec,wanc_bandpass)
xline(quakeSec)
xlabel("Time (seconds)")