write matlab code that will extract the peak-to-peak amplitude of the m-wave and h-reflex waveforms at each stimulation step with filtering and generating the variable in matlab

main.m
% Load EMG data
load emg_data.mat

% Apply a high-pass filter to remove baseline noise
high_pass_cutoff = 20; % Hz
sampling_rate = 1000; % Hz
[b, a] = butter(2, high_pass_cutoff / (sampling_rate / 2), 'high');
emg_data_filt = filtfilt(b, a, emg_data);

% Rectify the EMG signal
emg_data_rect = abs(emg_data_filt);

% Smooth the rectified EMG signal with a moving average filter
window_size = round(0.1 * sampling_rate); % 100 ms window
emg_data_smooth = movmean(emg_data_rect, window_size);

% Find the peaks and troughs of the smoothed EMG signal
[~, m_wave_pk_locs] = findpeaks(emg_data_smooth(1:stim_onset_loc), 'MinPeakDistance', round(0.01 * sampling_rate));
[~, h_reflex_pk_locs] = findpeaks(emg_data_smooth(stim_onset_loc:end), 'MinPeakDistance', round(0.01 * sampling_rate));
h_reflex_pk_locs = h_reflex_pk_locs + stim_onset_loc - 1;

[~, m_wave_tr_locs] = findpeaks(-emg_data_smooth(1:stim_onset_loc));
[~, h_reflex_tr_locs] = findpeaks(-emg_data_smooth(stim_onset_loc:end));
h_reflex_tr_locs = h_reflex_tr_locs + stim_onset_loc - 1;

% Calculate the peak-to-peak amplitude of each waveform
m_wave_amp = zeros(n_stim_steps, 1);
h_reflex_amp = zeros(n_stim_steps, 1);
for i = 1:n_stim_steps
    m_wave_pk = max(emg_data_filt(m_wave_pk_locs(i):m_wave_tr_locs(i)));
    m_wave_tr = min(emg_data_filt(m_wave_tr_locs(i):m_wave_pk_locs(i+1)));
    m_wave_amp(i) = abs(m_wave_pk - m_wave_tr);

    h_reflex_pk = max(emg_data_filt(h_reflex_pk_locs(i):h_reflex_tr_locs(i)));
    h_reflex_tr = min(emg_data_filt(h_reflex_tr_locs(i):h_reflex_pk_locs(i+1)));
    h_reflex_amp(i) = abs(h_reflex_pk - h_reflex_tr);
end

% Generate variable with peak-to-peak amplitude values for each waveform
peak_to_peak = table(m_wave_amp, h_reflex_amp);
1736 chars
41 lines