2009/3/30 João Luís Silva <jsi...@fc.up.pt>: > Hi, > > I wrote a script to calculate the *optical* autocorrelation of an > electric field. It's like the autocorrelation, but sums the fields > instead of multiplying them. I'm calculating > > I(tau) = integral( abs(E(t)+E(t-tau))**2,t=-inf..inf)
You may be in trouble if there's cancellation, but can't you just rewrite this as E(t)**2+E(t-tau)**2-2*E(t)*E(t-tau)? Then you have two O(n) integrals and one standard autocorrelation... Anne > with script appended at the end. It's too slow for my purposes (takes ~5 > seconds, and scales ~O(N**2)). numpy's correlate is fast enough, but > isn't what I need as it multiplies instead of add the fields. Could you > help me get this script to run faster (without having to write it in > another programming language) ? > > Thanks, > João Silva > > #-------------------------------------------------------- > > import numpy as np > #import matplotlib.pyplot as plt > > n = 2**12 > n_autocorr = 3*n-2 > > c = 3E2 > w0 = 2.0*np.pi*c/800.0 > t_max = 100.0 > t = np.linspace(-t_max/2.0,t_max/2.0,n) > > E = np.exp(-(t/10.0)**2)*np.exp(1j*w0*t) #Electric field > > dt = t[1]-t[0] > t_autocorr=np.linspace(-dt*n_autocorr/2.0,dt*n_autocorr/2.0,n_autocorr) > E1 = np.zeros(n_autocorr,dtype=E.dtype) > E2 = np.zeros(n_autocorr,dtype=E.dtype) > Ac = np.zeros(n_autocorr,dtype=np.float64) > > E2[n-1:n-1+n] = E[:] > > for i in range(2*n-2): > E1[:] = 0.0 > E1[i:i+n] = E[:] > > Ac[i] = np.sum(np.abs(E1+E2)**2) > > Ac *= dt > > #plt.plot(t_autocorr,Ac) > #plt.show() > > #-------------------------------------------------------- > > _______________________________________________ > Numpy-discussion mailing list > Numpy-discussion@scipy.org > http://mail.scipy.org/mailman/listinfo/numpy-discussion > _______________________________________________ Numpy-discussion mailing list Numpy-discussion@scipy.org http://mail.scipy.org/mailman/listinfo/numpy-discussion