FFT error 'not supported to carry out script fft as a function'

I want to plot a graph as below. so I wrote a program using fft. but error message 'not supported to carry out script fft as a function' displayed. What should I do?

syms t f
T=5.0*10^(-10);
roll = 0.3;%roll-off β
A = pi*t/T;
x(t)= sin(A)/A*cos(roll*A)/(1-(2*roll*t/T)^2);
ht=matlabFunction(x(t));
y=fft(x(t));

X = f;
Y = y;
plot(X,Y);

formula of x(t),X(f) and graph I want to plot(green line) are shown as follows

Answer 1

[Paul]

Hi 柊介 小山内,

fourier can return a closed form expression with a little help.

syms t w f real
T = sym(5.0)*10^(-10);
roll = sym(0.3);%roll-off β
A = sym(pi)*t/T;
x(t) = sin(A)/A*cos(roll*A)/(1-(2*roll*t/T)^2);

rewrite x(t) in terms of expoentials before taking the Fourier transform.

X(w) = simplify(fourier(rewrite(x(t),'exp'),t,w))

X(w) = 

Convert to Thz

syms fThz
X(fThz) = X(2*sym(pi)*(fThz*1e12))

X(fThz) = 

The plot doesn't look like yours, actually it looks like one cycle of yours. However, I also don't see how the the plots you've posted for X(f) match the equation you've posted for X(f)

xfunc = matlabFunction(X(fThz)/T);
figure
plot(-0.01:.00001:0.01,abs(xfunc(-0.01:.00001:0.01)))

[Walter Roberson]

Good use of rewrite()!

[Paul]

Thanks. But shouldn't have to. Don't understand why fourier() can't figure that out. There are also other identies that fourier() could have used to solve the problem.

[柊介 小山内]

Probably this graph can write to keep shfifting Paul's graph at 50GHz(0.05THz).

If I can define formula as below, I can write same graph as mentioned above.

y(fTHz+0.05)=y(fTHz)

but I don't know how to define this formula.

I wrote this folmula and error message 'invalid define of function or indexing. ' displayed.

What should I do?

[Paul]

Define everything fully parameterically and use the actual notation in the problem statement.

syms t w f T beta real
%T = sym(5.0)*10^(-10);
%beta = sym(0.3);%roll-off β
%A = sym(pi)*t/T;
%x(t) = sin(A)/A*cos(roll*A)/(1-(2*roll*t/T)^2);
Pi = sym(pi);
x(t) = sin(Pi*t/T)/(Pi*t/T)*cos(Pi*beta*t/T)/(1 - (2*beta*t/T)^2);
X(w) = simplify(fourier(rewrite(x(t),'exp'),t,w));
X(f) = X(2*sym(pi)*f);

Verify that imag(X(f)) = 0

simplify(imag(rewrite(X(f),'sincos')))

ans = 

0

So X(f) is just the real part:

X(f) = simplify(real(rewrite(X(f),'sincos')),100);

Compare the Matlab solution from fourier() with the given, closed-form solution.

Xclosed(f) = piecewise( ...
    0 <= abs(f) < (1-beta)/2/T, T, ...
    (1-beta)/2/T <= abs(f) < (1+beta)/2/T, T/2*(1 + cos(Pi*T/beta*(abs(f)-(1-beta)/2/T))), ...
    0);
figure;
fplot(subs([X(f) Xclosed(f)],[T beta],[1 0.3]),[-1 1])

Looks good.

Convert to THz and sub in the actual values for T and beta

syms fTHz
Tval = 5e-10;
X(fTHz) = subs(X(fTHz*1e12),[T beta],[Tval 0.3]);

Sum up three shifted copies of X(fTHz)

xfunc = matlabFunction(X(fTHz)/Tval);
fval = -0.01:.00001:0.01;
figure
plot(fval,(xfunc(fval) + xfunc(fval-0.003) + xfunc(fval+0.003)))

Or, more compactly

figure
plot(fval,(sum(xfunc(fval+[0 0.003 -0.003].'))))

I have no idea if this approach is going to yield the correct curve in a mathematical sense, because I have no idea what that green curve really is. At a minimum, I do see that the y-axis is labeled PSD, which suggests some sort of squaring operation. Hopefully this Answer has provided a framework for you to get what you need, whatever that may be.

[柊介 小山内]

I have define the graph as folmula and the other part of graph can write only shifting graph at 50GHz(0.05) intervals.

in above folmula, you can write center part of graph(six one from left) and The other part of graph can express as below.

right part of graph can express in this folmula and left part can express using '-n' instead of n.

n is integer and n is depend on how many mountain part the graph exist(In shown graph n=-5:5).

if you know, Please tell me how to storage each X(f) and plot X(f).

[Paul]

As shown above, define X(f), without the frequency offset, in a single function, call the function with offset frequency inputs, and then sum the results. At least that's how I'm interpreting what you want.

fval = -0.01:.00001:0.01;
figure
T = 5e-10;
beta = 0.3;
figure
plot(fval,(xfunc(fval,T,beta) + xfunc(fval-0.003,T,beta) + xfunc(fval+0.003,T,beta)))

function X = xfunc(fThz,T,beta)
X = 0*fThz;
f = fThz*1e12;
X(abs(f) <= (1-beta)/2/T) = T;
ii = (1-beta)/2/T < abs(f) & abs(f) < (1+beta)/2/T;
X(ii) = T/2*(1 + cos(pi*T/beta*(abs(f(ii))-(1-beta)/2/T)));
end

[柊介 小山内]

I have to integrate X(f) as follow.

in folmula(17),Gwdm is equivalent to X(f). I don't know about matlab well, so I don't know whether X(f) can integrate with Paul's definition and how can integrate X(f) if it can integrate. I'm sorry for my lack of explanation and spending your time for me. But my life depend on this question. So please help me.

[柊介 小山内]

I don't know whether this is useful for you, my program I wrote part of it published here.

syms t f fTHz w
T=1/(5.0*10^10);
roll = 0.3;%roll-off β
a=(1-roll)/(2*T)
A = pi*t/T;
%x(t)= sin(A)/A*cos(roll*A)/(1-(2*roll*t/T)^2);
nch = 11;%number of mountain part
m =T/2*(1+cos(pi*T/roll))*(f-a);
A=zeros(size(nch));
for n =  -round(nch/2):round(nch/2)
    f= -a+50*10^9*n:1.0*10^9:a+50*10^9*n　
   x(f)=piecewise(-a+50*10^9*n<=f<=a+50*10^9*n,T,a+50*10^9*n<=f<=(1+roll)/(2*T)+50*10^9*n,m);
   A(n)=x(f);
end

[柊介 小山内]

sorry I can solve my problem. thank you for help me!

Answer 2

[Walter Roberson]

you named your file fft.m which makes it impossible to call the Mathworks fft function. You need to rename your fft.m

[柊介 小山内]

I renamed fft.m. but in my program, error message 'invalid data type' displayed. I want to do fft without changing folmula about x(t). How I can fix this problem?

[Walter Roberson]

In your code, your x(t) is a symbolic function.

You can never fft() a symbolic function or symbolic expression. fft() is strictly the discrete fourier transform, which must be applied to numeric values (especially double precision), not to symbolic values (even if they are symbolic numbers).

If you wanted to take the fourier transform of a symbolic function or symbolic expression you would use fourier

[柊介 小山内]

I can carry out fourier transform. but I can't plot graph about fourier transform. Is it impossible to plot fourier transform?

syms t f
T=5.0*10^(-10);
roll = 0.3;%roll-off β
A = pi*t/T;
x(t)= sin(A)/A*cos(roll*A)/(1-(2*roll*t/T)^2);
y(f)=fourier(x(t));
ht = matlabFunction(y)

X = -0.5:0.01:0.5;
Y = ht(X*1.0*10^12);
plot(X,Y);
xlabel('f[THz]');
ylabel('PSD A.U');

in this program, error message 'function or parameter t don't recognized' displayed. I thought t was erased when fourier transform carried out, but didn't erased although doing fourier transform.

I want to plot X(f) and we can obtain X(f) by doing fourier transform x(t).

[柊介 小山内]

sorry I forgot to add sample. If graph plot carried out correctly, X(f) is ploted as below.(green line)

You don't worry about amplitude. Shape of graph is most important. This graph shows how much larger is the GNLI(blue line) than X(f).

[Walter Roberson]

syms t f
T = sym(5.0)*10^(-10);
roll = sym(0.3);%roll-off β
A = sym(pi)*t/T;
x(t)= sin(A)/A*cos(roll*A)/(1-(2*roll*t/T)^2);
y(f)=fourier(x(t))

y(f) = 

char(x)

ans = '-(sin(2000000000*t*pi)*cos(600000000*t*pi))/(2000000000*t*pi*(1440000000000000000*t^2 - 1))'

Notice that the result has unevaluated calls to fourier(). That means that fourier() was unable to compute the fourier transform of that function.

I checked on Wolfram Alpha, which was able to come up with a transformation... but MATLAB is not able to do so.

https://www.wolframalpha.com/input?i=fourier+-%28sin%282000000000*t*pi%29*cos%28600000000*t*pi%29%29%2F%282000000000*t*pi*%281440000000000000000*t%5E2+-+1%29%29