Using trapz to find integral of a function in polar coordinates

I have an xy grid over which I define two cosine functions and try to find the integral of their product using trapz. Given the orthogonality of cosine function, i.e., integral of cos(n*theta)*cos(m*theta) = pi if n = m and 0 if n is not equal to m. However, using the following code, the results do not agree with this property. Can someone suggest me a way of doing this?

if true
x = linspace(-1,1,1000);
y = linspace(-1,1,1000);
[X, Y] = meshgrid(x, y);
[theta, r] = cart2pol(X, Y); %cartesian to polar coordinates
step = theta(1,2)-theta(1,1);
n = 4;
m = 8;
E1 = cos(n*theta);
E2 = cos(m*theta);
a = step*trapz(step*trapz(E1.*E2));
end

In this case a is not zero.