假如你是一个资深的Matlab的专家，请用Matlab写一段提取PSCAD中电流ia的数据的FFT谐波分析程序

作为MATLAB专家，我可以为您提供一个基本的示例程序，用于从PSCAD仿真结果中提取电流ia的时间序列，并对其进行快速傅立叶变换（FFT）进行谐波分析。请注意，这假设您已经将PSCAD的仿真结果导出为CSV文件，其中包含电流数据。我们将使用`readmatrix`函数读取数据并使用`fft`函数进行分析。

首先，安装PSCAD MATLAB接口，如果尚未安装，可以参考PSCAD官方文档。然后，假设您的数据文件名为`current_data.csv`，并且电流列名为`ia`：

% 导入PSCAD仿真数据
data = readmatrix('current_data.csv');

% 确保我们只关注ia列
ia_data = data(:, find(strcmp(data.Properties.VariableNames, 'ia')));

% 提取时间戳
time_vector = data(:, 1);

% 对电流数据进行采样频率同步（如果需要）
% time_vector = linspace(0, size(ia_data, 1)/fs, size(ia_data, 1)); % fs是采样率

% 提取信号（如果时间序列有多个通道，可以选择一个进行FFT）
ia_signal = ia_data(1, :);

% 快速傅立叶变换
freq_domain_ia = fft(ia_signal);
freq = linspace(0, fs/2, length(freq_domain_ia)/2 + 1); % fs是采样频率

% 计算每个频率点的幅值和相位
amp = abs(freq_domain_ia(1:length(freq_domain_ia)/2+1));
phase_rad = angle(freq_domain_ia(1:length(freq_domain_ia)/2+1));

% 显示前几个主导谐波及其频率
harmonics = 1:5; % 调整此行以查看不同数量的谐波
for h = harmonics
    freq_index = floor(length(freq) * (h - 1) / (length(freq) - 1));
    fprintf('Harmonic %d: Frequency = %.2f Hz, Amplitude = %.2f\n', h, freq(freq_index), amp(freq_index));
end

% 如果需要绘制频谱图
figure;
plot(freq, amp);
xlabel('Frequency [Hz]');
ylabel('Amplitude');
title('Current IA FFT Spectrum');