Error in double integration

I need to evaluate a double integration but showing some error

lam=532*10^-9;
z=100;
k=2*pi/lam;
omega=30;s=5;
r=linspace(0,100,100);
w0=0.002;m=1;rho=1;p=1;
E = 1./(w0^2) + (1i*k)./(2*z);
Con1=(1i./(2*lam*z))*exp((omega^2./2 + r.^2 + omega)./E).*exp((-1i.*k.*r.^2)./(2*z)).*(pi./E).*(1./(2*1i*sqrt(E)))^m;
Con2=(1i./(2*lam*z))*exp((omega^2./2 + r.^2 + omega)./conj(E)).*exp((-1i.*k.*r.^2)./(2*z)).*(pi./conj(E)).*(1./(2*1i*sqrt(conj(E)))).^m;
syms ph phi

for l=0:m
    E0 = Con1.*((exp(-r.*(cos(ph)+sin(ph))).*hermiteH(l,1i*(omega./2 - r*cos(ph))./sqrt(E)).*hermiteH(m-l,1i*(omega./2 - r*sin(ph))./sqrt(E))) - (exp(r.*(cos(ph)+sin(ph))).*hermiteH(l,-1i*(omega./2 + r*cos(ph))./sqrt(E)).*hermiteH(m-l,-1i*(omega./2 + r*sin(ph))./sqrt(E))));
    E0s = Con2.*((exp(-r.*(cos(phi)+sin(phi))).*hermiteH(l,1i*(omega./2 - r*cos(phi))./sqrt(conj(E))).*hermiteH(m-l,1i*(omega./2 - r*sin(phi))./sqrt(conj(E)))) - (exp(r.*(cos(phi)+sin(phi))).*hermiteH(l,-1i*(omega./2 + r.*cos(phi))./sqrt(conj(E))).*hermiteH(m-l,-1i*(omega./2 + r*sin(phi))./sqrt(conj(E)))));
    I = ((1i*p).^(m-l)).*nchoosek(m,l).*E0.*E0s.*exp(-1i.*s.*(ph-phi)).*exp(((-2.*r.^2)-(2.*r.^2.*cos(phi-ph)))./(rho.^2));
end
IQ = I;

ph_min = 0;
ph_max = 2*pi;
phi_min = 0;
phi_max = 2*pi;

fun = matlabFunction(IQ,'Vars',[ph,phi]);
result = integral2(fun, ph_min, ph_max, phi_min, phi_max);

Error using integral2Calc>integral2t/tensor
Integrand output size does not match the input size.

Error in integral2Calc>integral2t (line 55)
[Qsub,esub] = tensor(thetaL,thetaR,phiB,phiT);

Error in integral2Calc (line 9)
[q,errbnd] = integral2t(fun,xmin,xmax,ymin,ymax,optionstruct);

Error in integral2 (line 105)
Q = integral2Calc(fun,xmin,xmax,yminfun,ymaxfun,opstruct);

Answer 1

[Torsten]

syms r ph phi
lam=532*10^-9;
z=100;
k=2*pi/lam;
omega=30;
w0=0.002;m=1;rho=1;p=1;
E = 1./(w0^2) + (1i*k)./(2*z);
Con1=(1i./(2*lam*z))*exp((omega^2./2 + r.^2 + omega)./E).*exp((-1i.*k.*r.^2)./(2*z)).*(pi./E).*(1./(2*1i*sqrt(E)))^m;
Con2=(1i./(2*lam*z))*exp((omega^2./2 + r.^2 + omega)./conj(E)).*exp((-1i.*k.*r.^2)./(2*z)).*(pi./conj(E)).*(1./(2*1i*sqrt(conj(E)))).^m;


for l=0:m
    E0 = Con1.*((exp(-r.*(cos(ph)+sin(ph))).*hermiteH(l,1i*(omega./2 - r*cos(ph))./sqrt(E)).*hermiteH(m-l,1i*(omega./2 - r*sin(ph))./sqrt(E))) - (exp(r.*(cos(ph)+sin(ph))).*hermiteH(l,-1i*(omega./2 + r*cos(ph))./sqrt(E)).*hermiteH(m-l,-1i*(omega./2 + r*sin(ph))./sqrt(E))));
    E0s = Con2.*((exp(-r.*(cos(phi)+sin(phi))).*hermiteH(l,1i*(omega./2 - r*cos(phi))./sqrt(conj(E))).*hermiteH(m-l,1i*(omega./2 - r*sin(phi))./sqrt(conj(E)))) - (exp(r.*(cos(phi)+sin(phi))).*hermiteH(l,-1i*(omega./2 + r.*cos(phi))./sqrt(conj(E))).*hermiteH(m-l,-1i*(omega./2 + r*sin(phi))./sqrt(conj(E)))));
    I = ((1i*p).^(m-l)).*nchoosek(m,l).*E0.*E0s.*exp(-1i.*l.*(ph-phi)).*exp(((-2.*r.^2)-(2.*r.^2.*cos(phi-ph)))./(rho.^2));
end
IQ = I;

ph_min = 0;
ph_max = 2*pi;
phi_min = 0;
phi_max = 2*pi;

fun = matlabFunction(IQ,'Vars',[ph,phi,r]);
r_array = linspace(0,10,100);
result = arrayfun(@(r)integral2(@(ph,phi)fun(ph,phi,r),ph_min, ph_max, phi_min, phi_max),r_array)

result = 
  Columns 1 through 10

  -0.0000 + 0.0000i   0.0000 + 0.0000i   0.0000 - 0.0000i   0.0000 - 0.0000i  -0.0000 + 0.0000i  -0.0000 + 0.0000i  -0.0000 - 0.0000i  -0.0000 + 0.0000i  -0.0000 - 0.0000i  -0.0000 + 0.0000i

  Columns 11 through 20

  -0.0000 + 0.0000i   0.0000 - 0.0000i  -0.0000 - 0.0000i   0.0000 + 0.0000i   0.0000 - 0.0000i   0.0000 + 0.0000i  -0.0000 - 0.0000i  -0.0000 - 0.0000i   0.0000 + 0.0000i   0.0000 + 0.0000i

  Columns 21 through 30

  -0.0000 + 0.0000i   0.0000 + 0.0000i  -0.0000 + 0.0000i   0.0000 - 0.0000i   0.0000 - 0.0000i  -0.0000 + 0.0000i  -0.0000 + 0.0000i  -0.0000 - 0.0000i  -0.0000 + 0.0000i  -0.0000 - 0.0000i

  Columns 31 through 40

   0.0000 + 0.0000i   0.0000 + 0.0000i  -0.0000 - 0.0000i  -0.0000 - 0.0000i   0.0000 - 0.0000i  -0.0000 - 0.0000i  -0.0000 - 0.0000i  -0.0000 + 0.0000i  -0.0000 + 0.0000i   0.0000 - 0.0000i

  Columns 41 through 50

  -0.0000 - 0.0000i   0.0000 - 0.0000i  -0.0000 - 0.0000i   0.0000 - 0.0000i  -0.0000 + 0.0000i  -0.0000 + 0.0000i   0.0000 - 0.0000i  -0.0000 - 0.0000i   0.0000 - 0.0000i  -0.0000 - 0.0000i

  Columns 51 through 60

   0.0000 - 0.0000i  -0.0000 + 0.0000i  -0.0001 + 0.0001i   0.0000 - 0.0001i  -0.0001 - 0.0001i   0.0002 - 0.0001i  -0.0002 - 0.0002i  -0.0000 - 0.0003i  -0.0001 + 0.0004i  -0.0001 + 0.0005i

  Columns 61 through 70

  -0.0003 - 0.0007i  -0.0009 - 0.0001i   0.0003 - 0.0012i  -0.0016 + 0.0002i  -0.0017 - 0.0012i   0.0017 + 0.0021i   0.0027 + 0.0023i  -0.0047 - 0.0003i  -0.0032 + 0.0053i  -0.0068 - 0.0043i

  Columns 71 through 80

  -0.0013 + 0.0105i  -0.0099 + 0.0098i   0.0154 - 0.0098i   0.0166 - 0.0174i  -0.0018 + 0.0315i   0.0370 + 0.0187i  -0.0221 + 0.0498i   0.0713 - 0.0064i   0.0819 + 0.0464i  -0.0977 - 0.0761i

  Columns 81 through 90

  -0.1530 - 0.0558i   0.1956 - 0.0875i  -0.0220 - 0.2810i   0.3640 - 0.0717i  -0.3396 - 0.3510i  -0.1472 - 0.6259i   0.1278 + 0.8370i   0.5641 + 0.9619i  -1.4387 - 0.2965i  -0.8911 + 1.7179i

  Columns 91 through 100

  -2.4397 - 0.7423i   1.1402 + 3.1615i  -0.5514 + 4.3955i   0.7996 - 5.7850i  -2.3376 - 7.3367i   9.5097 + 3.5597i   7.1866 -11.3004i  16.7884 + 5.4952i  -8.9062 -21.5352i   0.3631 -30.7457i