% =========================================================================
% GENERATE SAMPLE P300 DATA
% =========================================================================
% Instructor Script: Generate realistic P300 ERP data for student activity
%
% This script creates a synthetic but realistic P300 waveform that students
% can use for the ERP Peak Detective activity.
%
% Run this script to create: P300_data.mat
% =========================================================================

clear all
close all
clc

%% PARAMETERS
% -------------------------------------------------------------------------
% These control the characteristics of the synthetic ERP

% Time parameters
sampling_rate = 250; % Hz (250 samples per second)
epoch_start = -200;  % ms before stimulus
epoch_end = 800;     % ms after stimulus

% Component parameters
% P300 component
p300_latency = 350;   % ms (peak time)
p300_amplitude = 8;   % μV (peak height)
p300_width = 100;     % ms (width of the peak)

% N200 component (comes before P300)
n200_latency = 220;   % ms
n200_amplitude = -3;  % μV (negative)
n200_width = 80;      % ms

% P100 component (early visual response)
p100_latency = 100;   % ms
p100_amplitude = 2;   % μV
p100_width = 40;      % ms

% Noise level
noise_level = 0.5;    % μV (realistic background noise)

%% CREATE TIME VECTOR
% -------------------------------------------------------------------------
time = epoch_start:(1000/sampling_rate):epoch_end;

%% GENERATE COMPONENTS
% -------------------------------------------------------------------------
% Each component is modeled as a Gaussian function

% P300 component (main peak students will find)
p300 = p300_amplitude * exp(-((time - p300_latency).^2) / (2 * p300_width^2));

% N200 component (earlier negative peak)
n200 = n200_amplitude * exp(-((time - n200_latency).^2) / (2 * n200_width^2));

% P100 component (early positive deflection)
p100 = p100_amplitude * exp(-((time - p100_latency).^2) / (2 * p100_width^2));

% Combine all components
voltage = p300 + n200 + p100;

% Add realistic noise
voltage = voltage + noise_level * randn(size(voltage));

% Add a slight baseline drift (very common in real ERPs)
baseline_drift = 0.5 * sin(2*pi*time/1000) - 0.5;
voltage = voltage + baseline_drift;

%% ADD METADATA
% -------------------------------------------------------------------------
electrode_name = 'Pz';  % Typical electrode for P300
condition = 'Target';    % Oddball paradigm condition

%% SAVE DATA
% -------------------------------------------------------------------------
save('P300_data.mat', 'time', 'voltage', 'electrode_name', 'condition')

fprintf('=================================================\n')
fprintf('Sample P300 data created successfully!\n')
fprintf('=================================================\n')
fprintf('File saved: P300_data.mat\n')
fprintf('Variables included:\n')
fprintf('  - time: %d samples from %.0f to %.0f ms\n', length(time), time(1), time(end))
fprintf('  - voltage: ERP amplitude in μV\n')
fprintf('  - electrode_name: %s\n', electrode_name)
fprintf('  - condition: %s\n', condition)
fprintf('\nExpected peak measurements:\n')
fprintf('  Peak latency: ~%.0f ms\n', p300_latency)
fprintf('  Peak amplitude: ~%.1f μV\n', p300_amplitude)
fprintf('=================================================\n')

%% PLOT FOR VERIFICATION
% -------------------------------------------------------------------------
figure('Name', 'Generated Sample Data', 'Position', [100, 100, 900, 600])

subplot(2,1,1)
plot(time, voltage, 'b-', 'LineWidth', 1.5)
xlabel('Time (ms)')
ylabel('Amplitude (μV)')
title('Complete Synthetic P300 ERP')
grid on
xline(0, '--k', 'Stimulus', 'LineWidth', 1.5)
yline(0, '--k', 'Baseline', 'LineWidth', 1.5)

subplot(2,1,2)
plot(time, p300, 'r-', 'LineWidth', 2, 'DisplayName', 'P300')
hold on
plot(time, n200, 'b-', 'LineWidth', 2, 'DisplayName', 'N200')
plot(time, p100, 'g-', 'LineWidth', 2, 'DisplayName', 'P100')
xlabel('Time (ms)')
ylabel('Amplitude (μV)')
title('Individual Components')
legend('Location', 'best')
grid on
xline(0, '--k', 'Stimulus', 'LineWidth', 1.5)
yline(0, '--k', 'Baseline', 'LineWidth', 1.5)
hold off

%% OPTIONAL: CREATE ALTERNATIVE DATASETS
% -------------------------------------------------------------------------
% Uncomment this section to create multiple datasets with variations

% % Dataset 2: Smaller P300 (low amplitude)
% voltage2 = (p300 * 0.6) + n200 + p100 + noise_level * randn(size(voltage));
% save('P300_data_small.mat', 'time', 'voltage2', 'electrode_name', 'condition')
% 
% % Dataset 3: Later P300 (slower processing)
% p300_late = p300_amplitude * exp(-((time - 450).^2) / (2 * p300_width^2));
% voltage3 = p300_late + n200 + p100 + noise_level * randn(size(voltage));
% save('P300_data_late.mat', 'time', 'voltage3', 'electrode_name', 'condition')
% 
% % Dataset 4: Noisier data (more realistic)
% voltage4 = p300 + n200 + p100 + (noise_level*2) * randn(size(voltage));
% save('P300_data_noisy.mat', 'time', 'voltage4', 'electrode_name', 'condition')

fprintf('\nData generation complete! Ready for student activity.\n')
