信号求和，同时保留每个单独时刻的相位矩阵的随机值

Question

当

t = const

时，我将频率范围从 0.8 到 11.1 的信号的所有谐波相加，步长为 0.1。每个单独的谐波（在范围内具有其自己单独频率的每个信号）都有其自己的随机相位。这是程序：

clear,clc
%% Data 
x = 0:1:100; 
a = 5.76; 
omega = 0.8:0.1:11.1; 
k = omega / a; 
phase = 2*pi*rand(length(omega),1);
for i = 1:length(x)
    PHASE(i,:) = phase(:,1);
end

[X,OMEGA] = meshgrid(x,omega);
K = OMEGA / a; 

t = 0; 

y = cos(OMEGA * t - K .* X + PHASE');

y_sum = sum(y);


figure(1);
plot(x, y(1:25:end,:));
title('y(x)');
grid on;

figure(2);
plot(x, y_sum);
title('y(x)');
grid on;

为每个单独的频率生成一个随机相位，但相对于

x

坐标，它不应该改变，这就是为什么我在开始时构建了这个循环，它构建了与这些条件相对应的随机相位矩阵（是否有任何不使用循环的其他方法可以做到这一点？）。

任务如下：对于每个单独的时刻

t

（现在有几个），对信号进行求和，但同时，以便存储

PHASE

矩阵的值对于每个后续时刻，并且不会以新的方式更新。

起初我尝试通过循环来做到这一点，但后来我意识到每次迭代

t

它都会更新

PHASE

矩阵的元素值。结果，图表不正确。

结果，我应该得到这样的图表： 只有我手动为每个单独的时刻分别注册函数

y

和

y_sum

，结果才会这样，这是非常不方便的（特别是如果我想要更改时间范围）。

t0 = 0; t1 = 0.1; t2 = 0.5; t3 = 1;

y0 = cos(OMEGA * t0 - K .* X + PHASE');
y1 = cos(OMEGA * t1 - K .* X + PHASE');
y2 = cos(OMEGA * t2 - K .* X + PHASE');
y3 = cos(OMEGA * t3 - K .* X + PHASE');

y_sum0 = sum(y0);
y_sum1 = sum(y1);
y_sum2 = sum(y2);
y_sum3 = sum(y3);

figure(2);
plot(x, [y_sum0',y_sum1',y_sum2',y_sum3']);
legend('t0 = 0 s', 't1 = 0.1 s', 't2 = 0.5 s', 't3 = 1 s', 'Location','eastoutside');
title('y(x)');
grid on;

如何以更智能的方式执行这样的任务，而不需要在各个时刻增加功能

y

和

y_sum

？

Answer 1

希望下面的示例代码能够回答您的问题。

几点建议：

1.) 就我个人而言，我从不从脚本中调用clear 或clc。虽然如果脚本在您自己的计算机上运行，使用这些命令非常有趣和游戏，但一旦您将脚本交给同事，您就不知道这些命令可能会删除同事没有删除的内容想要删除。

2.) 不要害怕将代码编写为函数而不是脚本。这样做可以划分代码，以便更容易理解。它还可以让您避免剩余变量使工作区变得混乱（这反过来会减少调用clear和clc的需要）。最后，如果您愿意，可以将函数保存在单独的文件中，以便将来可以重复使用它们，而无需复制和粘贴代码。

3.) 有时，循环确实是最好的做事方式，特别是当您发现自己正在处理两个以上维度的数据时。别想太多。在这种情况下，下面的示例代码可能远没有达到导致明显性能滞后的复杂性，即使它使用循环。

4.) 您的绘图越先进，您就会发现它对于保留您创建的图形对象的句柄就越有帮助。这样做可以让您在绘制数据后执行更改颜色、创建自定义图例标签等操作。此外，如果您想编写 MATLAB 应用程序，这几乎是处理应用程序内部绘图的唯一方法。

我在Notepad++中编写了下面的示例代码，因为我现在使用的计算机上没有安装MATLAB。不幸的是，这意味着我几乎可以肯定它将有一两个愚蠢的语法错误阻止它运行。不过，它应该非常接近。

% define the inputs for the problem
t = [0, 0.1, 0.5, 1];
omega = 0.8:0.1:11.1; 
x = 0:1:100;
a = 5.76;

% perform the calculation - the calcY function is declared at the bottom of the script
[y, y_sum] = calcY(t, omega, x, a);


% initialize the figure windows
yFigHandle = figure();                  % create a new figure
yAxesHandle = axes(yFigHandle);         % create a new set of cartesian axes within the new figure
title(yAxesHandle,'y(x)');
yLineObjects = gobjects(numel(y),1);    % preallocate a place holder array of graphics objects where we will eventually save the plotted lines
hold(yAxesHandle,'on');

y_sumFigHandle = figure();
y_sumAxesHandle = axes(y_sumFigHandle);
title(y_sumAxesHandle,'y_sum(x)');
y_sumLineObjects = gobjects(numel(y_sum),1);
hold(y_sumAxesHandle,'on');

% plot the data
for i = 1:numel(y)
    % By the way, I'm not sure why we're only plotting every 25th element here...there should be 114 rows in each
    % element of the yCell and that number is not divisible by 25. Once again, that's only if my back-of-the-napkin
    % math is correct, since I haven't been able to run this code.
    yLineObjects(i) = plot(yAxesHandle, x, y{i}(1:25:end,:));
    y_sumLineObjects(i) = plot(y_sumAxesHandle, x, y_sum{i});
end

% finish configuring the plots
hold(yAxesHandle,'off');
grid(yAxesHandle,'on');

hold(y_sumAxesHandle,'off');
grid(y_sumAxesHandle,'on');



function [yCell, y_sumCell] = calcY(tVector, omegaVector, xVector, a)
    % [yCell, y_sumCell] = calcY(tVector, omegaVector, xVector, a)
    %
    % Performs your calculation (which I'm still not entirely certain the purpose of, but that's not really important).
    % 
    % Returns two cell arrays, yCell and y_sumCell. They each have one column, and the same number of rows as elements
    % in tVector.
    %   Each element of yCell is a matrix of doubles with number of rows = number of elements in omegaVector and number
    %      of columns = number of elements in xVector.
    %   Each element of y_sumCell is a row vector of doubles with number of elements = number of elements in xVector.


    % ALWAYS validate your function inputs. This makes it clear to anyone else looking at this code what is expected by
    % this function, and it also prevents stupid bugs where you fed the function the wrong thing.
    arguments (Input)
        tVector (1,:) double {mustBeNonempty, mustBeReal, mustBeNonnegative}
        omegaVector (1,:) double {mustBeNonempty, mustBeReal, mustBeNonnegative}
        xVector (1,:) {mustBeNonempty, mustBeReal, mustBeNonnegative}
        a (1,1) double {mustBeReal}
    end

    % Ouput validation is optional, but I find it useful for the same reasons as above.
    arguments (Output)
        yCell (:,1) cell {mustBeNonempty}
        y_sumCell (:,1) cell {mustBeNonempty}
    end


    % determine the sizes of the output arrays
    numRows = numel(omegaVector);
    numColumns = numel(xVector);
    numCells = numel(tVector);

    % configure the inputs for the calculation
    omegaMatrix = omegaVector' .* ones(1,numColumns);       % make omegaVector a column vector and copy it to every column of omegaMatrix
    k = (omegaVector/a)';
    phase = 2*pi*rand(1,numRows);

    % preallocate the cell arrays that will be returned
    yCell = cell(numCells,1);
    y_sumCell = y;

    % loop through each element of tVector and perform the calculation.
    for i = 1:numCells
        % NOTE: If you perform element-wise multiplication between an N-element column vector and an M-element row vector,
        %   the result is an NXM matrix. Similarly, I believe that if you add a row vector to a matrix that has the same
        %   number of columns, the row vector will get added to every row.
        %   The subtraction operation is the one I'm not sure about, hence why I'm using omegaMatrix instead of just
        %   omegaVector'. Once again, I don't have MATLAB installed to test this, so I'm kind of guessing.
        yCell{i} = cos(omegaMatrix * tVector(i) - k .* xVector + phase);
        y_sumCell{i} = sum(yCell{i});
    end
end