写一段matlab代码，提取语音信号的mfcc特征（语音信号存储在audio_data中）

% 假设语音信号已经存储在audio_data中

% 设置参数
frame_length = 20; % 帧长为20ms
frame_shift = 10; % 帧移为10ms
num_mfcc = 12; % 提取12个MFCC系数
num_filters = 26; % 设置26个三角滤波器

% 将语音信号分帧
frame_size = frame_length*0.001*Fs; % 计算帧长的样点数
frame_shift_size = frame_shift*0.001*Fs; % 计算帧移的样点数
num_frames = floor((length(audio_data)-frame_size)/frame_shift_size)+1; % 计算总帧数
frames = zeros(frame_size,num_frames); % 存储所有帧
for i = 1:num_frames
start_index = (i-1)*frame_shift_size+1;
end_index = start_index+frame_size-1;
frames(:,i) = audio_data(start_index:end_index);
end

% 提取每帧的MFCC系数
mfcc = zeros(num_mfcc,num_frames); % 存储所有帧的MFCC系数
for i = 1:num_frames
% 对每帧进行预加重
pre_emphasis_factor = 0.97; % 预加重系数
pre_emphasis_frame = frames(:,i) - [0; frames(1:end-1,i)]*pre_emphasis_factor;
% 对每帧进行加窗
window = hamming(frame_size);
windowed_frame = pre_emphasis_frame.*window;
% 对每帧进行FFT
fft_size = 2^nextpow2(frame_size); % 选择FFT点数为2的幂次方
fft_frame = fft(windowed_frame,fft_size);
% 计算每帧的功率谱
power_spectrum = abs(fft_frame).^2/fft_size;
% 计算三角滤波器组
filter_bank = zeros(num_filters,fft_size/2+1); % 存储所有滤波器
% 设置三角滤波器组的三个顶点频率
low_freq = 0;
high_freq = Fs/2;
mel_low = 2595*log10(1+low_freq/700);
mel_high = 2595*log10(1+high_freq/700);
mel_points = linspace(mel_low,mel_high,num_filters+2);
hz_points = 700*(10.^(mel_points/2595)-1);
% 计算每个滤波器的通带
for j = 1:num_filters
filter_start = hz_points(j);
filter_peak = hz_points(j+1);
filter_end = hz_points(j+2);
filter_bank(j,:) = tri_filter(filter_start,filter_peak,filter_end,fft_size,Fs);
end
% 对每帧的功率谱进行滤波器组加权
filtered_spectrum = filter_bank*log10(power_spectrum);
% 对每帧的滤波器组加权功率谱取DCT系数
num_ceps = num_mfcc+1; % DCT系数个数
ceps = dct(filtered_spectrum);
mfcc(:,i) = ceps(2:num_ceps); % 不考虑第一个DCT系数
end

% 定义三角滤波器函数
function [filter_response] = tri_filter(filter_start,filter_peak,filter_end,fft_size,Fs)
filter_start_index = floor(filter_start/Hz*fft_size)+1;
filter_peak_index = floor(filter_peak/Hz*fft_size)+1;
filter_end_index = floor(filter_end/Hz*fft_size)+1;
filter_response = zeros(1,fft_size/2+1);
% 上升部分
for i = filter_start_index:filter_peak_index
filter_response(i) = (i-filter_start_index)/(filter_peak_index-filter_start_index);
end
% 下降部分
for i = filter_peak_index+1:filter_end_index
filter_response(i) = (filter_end_index-i)/(filter_end_index-filter_peak_index);
end
end