memory management

[Sudheer Joseph]

HI,
        I have been trying to compute cross correlation between a time series at a location f(1) and the timeseries of spatial data f(XYT) and saving the resulting correlation coefficients and lags in a 3 dimensional array which is of fairly big size. Though the code I made for this purpose works up to few iterations then it hangs due to apparent memory crunch. Can anybody suggest a better way to handle this situation so that the computation and data storing can be done with out hangups. Finally I intend to save the data as netcdf file which is not implemented as of now. Below is the piece of code I wrote for this purpose.

from mpl_toolkits.basemap import Basemap as bm, shiftgrid, cm
import numpy as np
import matplotlib.pyplot as plt
from netCDF4 import Dataset
from math import pow, sqrt
import sys
from scipy.stats import t
indep=120
nlags=365
ncin = Dataset('qu_ru.nc', 'r')
lons = ncin.variables['LON421_600'][:]
lats = ncin.variables['LAT81_220'][:]
dep = ncin.variables['DEPTH1_29'][:]
adep=(dep==indep).nonzero()
didx=int(adep[0])
qu = ncin.variables['qu'][:,:,:]
#qv = ncin.variables['QV'][0,:,:]
ru = ncin.variables['ru'][:,didx,0,0]
ncin.close()
fig = plt.figure()
ax = fig.add_axes([0.1,0.1,0.8,0.8])
# use major and minor sphere radii from WGS84 ellipsoid.
m = bm(projection='cyl', llcrnrlon=30, llcrnrlat=-40,urcrnrlon=120, urcrnrlat=30)
# transform to nx x ny regularly spaced 5km native projection grid
nx = int((m.xmax-m.xmin))+1; ny = int((m.ymax-m.ymin)+1)
q=ru[1:2190]
qmean=np.mean(q)
qstd=np.std(q)
qnorm=(q-qmean)/qstd
lags3d=np.arange(731*140*180).reshape(731,140,180)
r3d=np.arange(731*140*180).reshape(731,140,180)
for i in np.arange(len(lons)):
   for j in np.arange(len(lats)):
      print i,j
      p=qu[1:2190,j,i].squeeze()
      p.shape
      pmean=np.mean(p)
      pstd=np.std(p)
      pnorm=(p-pmean)/pstd
      n=len(p)
#      fg=plt.figure()
      c=plt.xcorr(p,q,usevlines=True,maxlags=nlags,normed=True,lw=2)
      acp=plt.acorr(p,usevlines=True,maxlags=nlags,normed=True,lw=2)
      acq=plt.acorr(q,usevlines=True,maxlags=nlags,normed=True,lw=2)
      acp[1][nlags]=0
      acq[1][nlags]=0
      lags=c[0]
      r=c[1]
      lags3d[:,j,i]=lags
      r3d[:,j,i]=r