Evaluating a complicated integral in Matlab gives NaN

I wish to plot the values of the following integral as a function of ψ.

where, and

Following is my code:

Nu = 4;             % Group cardinality
Pt = 10^(-12);      % Effective transmission power
fc = 3*10^9;        % Carrier frequency
df = 10^4;          % Sub-carrier spacing
c = 3*10^8;         % Speed of light
dv = 0.01;
Vmax = 0:dv:100-dv; % Maximum velocity of UEs
b = pi*Vmax*fc/c/df;

f2 = zeros(1, length(b));
f0 = @(t) sin(t)/t;    % Defining Si(z) function
for i = 1:length(b)    
    f1 = @(psi) integral(f0, 0, 2*b(i)*cos(psi), 'ArrayValued', true)./b(i)./cos(psi) - (sin(b(i)*cos(psi))).^2./b(i)^2./(cos(psi)).^2;
    f2(i) = integral(f1, 0, pi/2, 'ArrayValued', true);
    P_ICI = Nu*Pt - (2/pi)*Nu*Pt*f2;
end

I use the integral function of MATLAB to compute the integral but it gives me the following error.

Warning: Infinite or Not-a-Number value encountered.
In integralCalc/midpArea (line 397)
  In integralCalc (line 65)
  In integral (line 88)
  In @(psi)integral(f0,0,2*b(i)*cos(psi),'ArrayValued',true)./b(i)./cos(psi)-(sin(b(i)*cos(psi))).^2./b(i)^2./(cos(psi)).^2
  In integralCalc/iterateArrayValued (line 193)
  In integralCalc/vadapt (line 130)
  In integralCalc (line 75)
  In integral (line 88)
  In CA_1 (line 13) 

Can anyone suggest a way to get around this? Or point out what is wrong? Many thanks!

[darova]

Try small number instead of 0

[MaxFrost]

Thanks! I got it.

Answer 1

[Matt J]

Your integrand has a singularity at psi = pi/2. So perhaps try,

f2(i) = integral(f1, 0, pi/2*(1-eps), 'ArrayValued', true);