当前位置：服务支持 > 软件文章 > MATLAB编写m函数理解：y=f(g(x))*h(x)复合函数编程实现

MATLAB编写m函数理解：y=f(g(x))*h(x)复合函数编程实现

阅读数 2

任务：

将数学形式： $y=f(g(x))$ ，编写为MATLAB m 函数。

注：

1.任何字母本质上自变量都是k

MATLAB中，任何字母本质上自变量都是k，k表示第k个数。即任何函数再复杂，再复合，都可以for循环笨方法一个一个算出来。

$y=[y(k)]=[f(g( x(k) ))]$ ，y,x,[]表示矩阵，y(k),x(k)表示一个数。

则MATLAB里，

k=1:K

x=[x(k)]=x1:xs:x2

$y(k)=f(g(x(k)))$

+for 循环

end

但下标自变量k只能从1开始，不能从0。但字母值可以从零开始，如n=0:10;

2.除了笨方法，可能可以用特别的技巧：

2.1 简单函数作用矩阵=矩阵，代替for循环

若数学形式：y=f(x)，f都是简单的函数，比如exp,sin等，则MATLAB里可：

y=f(x) % y,x均为矩阵

不用

for k..

y(k)=f(x(k))

end

2.2 矩阵点乘替代函数相乘

若数学形式：y=f(x)g(x)，f,g都是简单的函数，比如exp,sin等，则MATLAB里或许可以省去for循环，通过：

y=f(x) .* g(x) % y,x均为矩阵。

2.3 矩阵相乘替代求和符号

若数学形式： $y=\sum_{k=1}^{N}f(x(k)))\cdot g((x(k)))$ ，则MATLAB里或许可以矩阵相乘替代求和符号，通过：

f=[f(k)]=[f(x(k))]

g=[g(k)]=[g(x(k))]

y=f * g'，f和g是k所控制的最高层变量。

3.例题：DTFT *Ts

求离散时间傅里叶变换 X(exp(jw))* Ts

编程出： $X(e^{j\Omega T_{s}})=\sum_{n=0}^{N-1}x(n)e^{-j\Omega T_{s}n} \cdot T_{s}$ 。

解：

1. X自变量看似exp(j*omega*T)，很复杂，但其实是X=[X(k)]。

采用逻辑自底层而上。

写成好看的形式：

Ω(𝑘)=Ω1:Δ:Ω2 Ω( k )=Ω1:Δ:Ω2 \Omega ({\color{Red} k})=\Omega1:\Delta :\Omega2

𝑥=[𝑥(𝑘)]=𝑒−𝑗Ω𝑇𝑠𝑛(𝑘)𝑐𝑜𝑠(𝑤0∗𝑇∗𝑛(𝑘))) x=[x( k )]=e − j Ω T s n ( k ) cos(w0∗T∗n( k ))) x=[x({\color{Red} k})]=e^{-j\Omega T_{s} n({\color{Red} k})}cos(w0 * T*n({\color{Red} k})))。

𝑋(𝑘)=∑𝑁−1𝑛=0𝑥(𝑛(𝑘′)))𝑒−𝑗Ω(𝑘)𝑇𝑠𝑛(𝑘′)⋅𝑇𝑠 X( k )= ∑ n = 0 N − 1 x(n( k ′ )))e − j Ω ( k ) T s n ( k ′ ) ⋅ T s X({\color{Red} k})=\sum_{n=0}^{N-1}x(n({\color{Blue} k'})))e^{-j\Omega({\color{Red} k}) T_{s}n({\color{Blue} k'})} \cdot T_{s}其实k'没必要标了，只要写到自变量（k），即

𝑋(𝑘)=∑𝑁−1𝑛=0𝑥([𝑛])𝑒−𝑗Ω(𝑘)𝑇𝑠[𝑛]⋅𝑇𝑠 X( k )= ∑ n = 0 N − 1 x([n])e − j Ω ( k ) T s [ n ] ⋅ T s X({\color{Red} k})=\sum_{n=0}^{N-1}x([n])e^{-j\Omega({\color{Red} k}) T_{s}[n]} \cdot T_{s},

2. exp, cos 为简单函数，MATLAB里可以省去for，直接作用矩阵产生矩阵。

3. exp * cos 简单函数相乘，MATLAB里可以点乘 .* 。

𝑥=[𝑥(𝑘)]=[𝑒−𝑗Ω𝑇𝑠𝑛(𝑘)𝑐𝑜𝑠(𝑤0∗𝑇∗𝑛(𝑘)))]=𝑒−𝑗Ω𝑇𝑠[𝑛].∗𝑐𝑜𝑠(𝑤0∗𝑇∗[𝑛])) x = [ x ( k ) ] = [ e − j Ω T s n ( k ) c o s ( w 0 ∗ T ∗ n ( k ) ) ) ] = e − j Ω T s [ n ] . ∗ c o s ( w 0 ∗ T ∗ [ n ] ) ) \begin{array}{c} x = [x(k)] =[{e^{ - j\Omega {T_s}n({\color{Red} k})}}cos(w0*T*n({\color{Red} k})))]\\ {\rm{ = }}{e^{ - j\Omega {T_s}{\color{Green} [n]}}}.*cos(w0*T*{\color{Green} [n]})) \end{array}\:

4. 有求和符号，MATLAB可以通过矩阵相乘简化。

𝑋(𝑘)=∑𝑛=0𝑁−1𝑥(𝑛(𝑘′)))𝑒−𝑗Ω(𝑘)𝑇𝑠𝑛(𝑘′)⋅𝑇𝑠=𝑇𝑠𝑒−𝑗Ω(𝑘)𝑇𝑠[𝑛]∗[𝑥]′ X ( k ) = ∑ n = 0 N − 1 x ( n ( k ′ ) ) ) e − j Ω ( k ) T s n ( k ′ ) ⋅ T s = T s e − j Ω ( k ) T s [ n ] ∗ [ x ]′ \begin{array}{c} X({\color{Red} k}) = \sum\limits_{n = 0}^{N - 1} x (n({\color{Blue} k'}))){e^{ - j\Omega ({\color{Red} k}){T_s}n({\color{Blue} k'})}} \cdot {T_s}\\ {\rm{ = }}{T_s}{e^{ - j\Omega ({\color{Red} k}){T_s}[n]}}*[x]' \end{array}

代码：

% 复合函数% 求DTFT*Tsclc;clear;close all;  %% x% 参数fs=3;   %采样频率f0=0.5;omega0=2*pi*f0;Ts=1/fs;T=100; %采样时间N=T*fs;%自变量n=0:N-1;zeta=0.05;x=exp(-zeta*Ts*[n]) .* cos(omega0*Ts*[n]); %% Xomega=[0: 0.01: 4*pi] /Ts;    %omega不必与[x(k)]数目相同，因为DTFD的w是连续的for k=1:length(omega)     X(k)=Ts * exp(-1i*omega(k)*Ts*[n]) * [x]';end% 连续时间傅里叶变换 FTY=1/2*(1./(zeta+1i*(omega- omega0))+1./(zeta+1i*(omega+omega0))); h=figure(1);set(h,'Position',[150,150,800,500])semilogy(omega/2/pi,abs(Y),'b','LineWidth',2)xlabel('\Omega/2 \pi =f','FontName','Times New /*不要回车换行*/ Roman','FontSize',14)ylabel('Amplitude','FontName','Times New /*不要回车换行*/ Roman','FontSize',14)grid onhold on semilogy(omega/2/pi,abs(X),'k','LineWidth',2)set(gca,'FontName','Times New Roman','FontSize',14)legend('True,FT',['DTFT*Ts,f_s=',num2str(fs)])

仿真结果：