0

我对八度有一个小问题。

我想模拟一些东西，因此我需要一个循环，但不幸的是我的数据没有按八度保存。我尝试了一些解决这个问题的可能性，但我找不到任何解决方案。

这是我的文件的代码：

%------ input from the user --------
%-----------------------------------

at=1;                  % acq. time in s

nu=300;               % Resonance frequency in Hz

scantime=1.1;          % time for a scan, >at in s

T1=(0.1*at:at/10:120*at)';                 % in s

T2=T1;                      % in s

noisestd=1;                 % noise std

expt = (at+scantime);                % experiment time



%-------RATIOS----------


T2T1ratio=(0.01:0.01:12)';                        %T2T1 ratio


att1ratio=at*1./(T1);                                %acquisition time by T1


%----------CALL CALCSPEC.M--------------
%----------------------------------------



for T1idx = (1:(rows(T1)))'
    for T2T1ratioidx = (1:rows(T2T1ratio))'
        for att1ratioidx = (1:rows(att1ratio))'
            spectrum = [s2n, peakwidth] = calcspec(at,nu,scantime,T1,T2,noisestd,expt);
            m={T1idx,T2T1ratioidx,att1ratioidx,s2n,peakwidth}
        endfor
    endfor
endfor

这是我的 .m 文件的代码：

function [s2n, peakwidth] = calcspec(at,nu,scantime,T1,T2,noisestd,expt);



%---- calculating input vars -----------------
%---------------------------------------------

ns=round(expt/scantime);                % nr. of scans

Fs = 4*nu;                              % Sampling frequency (>=2*nu by Nyquist)

t = ((0:(at*Fs-1))/Fs)';                   % Time vector

fn=2^nextpow2(at*Fs);                   % number of points for the FT (must be power of 2)

res = Fs/fn;                            % spectral resolution (bins in the fft)

freq=0:res:(Fs-res);                    % frequency axis

fid=exp(-i*2*pi*nu.*t).*exp(-t./T2);    % generation of "clean" FID

fid=ns*fid'/max(fid);                   % normalization of FID 

for scan=1:ns
    noisyfid = fid+=noisestd*randn(1,res+1)';       % noisy fid
endfor




% FOURIER TRANSFORM 
% -----------------

sig = fft(fid, fn);                                                   % FT of clean spec

noisysig = fft(noisyfid, fn);                                        % Fourier Transform of noisy spec


%-------Calc S2N ----------
%---------------------

s2n = var(real(sig))/var(real(noisysig));


% DETERMINE FULL WIDTH AT HALF MAXIMUM + INTENSITY
%----------------------------------------------------

peakwidth = fwhm(real(sig));

intensity = max(real(sig));

如您所见（希望如此），我尝试将数据保存在单元格数组中，但计算后只保存了 s2n 和 peakwidth 的一个值。

有人可以告诉我我做错了什么吗？

4

3 回答 3

0

非常感谢您的回答！

我做的两个=标志是因为否则我的八度音程给我的信息是有一些变量是未定义的，所以我尝试了这个。

现在我试图改变所有错误的东西，但我仍然不明白索引的那个东西。现在我尝试T1idx在我的代码中使用索引，以便我现在得到：

    for T1idx = (1:(rows(T1)))
    for T2T1ratioidx = (1:rows(T2T1ratio))
        for att1ratioidx = (1:rows(att1ratio))
            [sig(T1idx),noisysig(T1idx)] = calcspec(at,nu,scantime,T1,T2,noisestd,expt);
            s2n = sig/noisysig;
            peakwidth = fwhm(real(sig));
            m={s2n,peakwidth};
        endfor
    endfor
endfor

但现在我得到了消息：A(I) = X: X must have the same size as I我没有任何好主意来获得“更好”的指数......

最后，我只想得到一个具有以下结构的单元格数组：

T1idx   T2T1ratioidx   att1ratioidx       s2n         peakwidth
1          1                1          a value         a value
2          .                .              .              .
3          .                .              .              .
4
.

很清楚我的意思是什么？

于 2013-06-24T13:49:47.997 回答
0

我认为以前一切都很好，除了换位。为了保存所有你可能做的事情

for T1idx = (1:(rows(T1)))
    for T2T1ratioidx = (1:rows(T2T1ratio))
        for att1ratioidx = (1:rows(att1ratio))
            spectrum = [s2n, peakwidth] = calcspec(at,nu,scantime,T1,T2,noisestd,expt);
            m=[m; {T1idx,T2T1ratioidx,att1ratioidx,s2n,peakwidth}]
        endfor
    endfor
endfor

编辑：让 m 看起来像你想要的那样：

i = 1;
m = [];
for T1idx = (1:(rows(T1)))
    for T2T1ratioidx = (1:rows(T2T1ratio))
        for att1ratioidx = (1:rows(att1ratio))
            spectrum = [s2n, peakwidth] = calcspec(at,nu,scantime,T1,T2,noisestd,expt);
            m(i,:)=[T1idx,T2T1ratioidx,att1ratioidx,s2n,peakwidth];
            i = i + 1;
        endfor
    endfor
endfor
m
于 2013-06-24T15:38:36.523 回答
0

尽管您正在循环，但您在循环内的计算中根本没有使用索引。你可能需要做这样的事情：

[s2n(%appropriate index%), peakwidth(%appropriate index%)] = calcspec(at,nu,scantime,T1(T2T1ratioidx),T2(T2T1ratioidx),noisestd,expt);
% same applies for m calculation

此外，您在同一行上有两个=标志，不确定最终结果是什么，但可能不是您所期望的：

spectrum = [s2n, peakwidth] = ...
于 2013-06-24T09:01:49.977 回答