Index: trunk/tools/nc_var_tools.py
===================================================================
--- trunk/tools/nc_var_tools.py	(revision 282)
+++ trunk/tools/nc_var_tools.py	(revision 283)
@@ -10404,5 +10404,6 @@
             trajvals[it,2] = latv[gtrajvals[it,2],gtrajvals[it,1]]
 
-#            print iline, it,'time:',trajvals[it,0],'lon lat:', trajvals[it,1], trajvals[it,2]
+#            print iline, it,'time:',trajvals[it,0],'lon lat:', trajvals[it,1],       \
+#              trajvals[it,2]
 
 # Assuming time as the first dimension in a fortran like way (t=0 as 1)
@@ -10411,5 +10412,5 @@
 
 # box values
-            if gtrajvals[it,2]-box2 < 0 or gtrajvals[it,2]+box2 + 1 > dimy + 1 \
+            if gtrajvals[it,2]-box2 < 0 or gtrajvals[it,2]+box2 + 1 > dimy + 1       \
               or gtrajvals[it,1]-box2 < 0 or gtrajvals[it,1]+box2 + 1 > dimx + 1:
 
@@ -10495,6 +10496,6 @@
                 cyrangeslice.append([gtrajvals[it,2]-Nrad, gtrajvals[it,2]+Nrad + 1])
                 cxrangeslice.append([gtrajvals[it,1]-Nrad, gtrajvals[it,1]+Nrad + 1])
-                cyrangeslice2D.append([gtrajvals[it,2]-Nrad, gtrajvals[it,2]+Nrad + 1])
-                cxrangeslice2D.append([gtrajvals[it,1]-Nrad, gtrajvals[it,1]+Nrad + 1])
+                cyrangeslice2D.append([gtrajvals[it,2]-Nrad, gtrajvals[it,2]+Nrad+1])
+                cxrangeslice2D.append([gtrajvals[it,1]-Nrad, gtrajvals[it,1]+Nrad+1])
             it = it + 1
             iline = iline + 1
@@ -10539,57 +10540,42 @@
             dimz = varobj.shape[1]
  
-            if not objofile.dimensions.has_key('z'):
-                newdim = objofile.createDimension('z', dimz)
-
-                varzobj = objfile.variables[zn]
-                newvar = objofile.createVariable(zn, 'f8', ('z'))
-
-                if searchInlist(varzobj.ncattrs(),'standard_name'):
-                    vsname = varzobj.getncattr('standard_name')
-                else:
-                    vsname = variables_values(zn)[1]
-                if searchInlist(varzobj.ncattrs(),'long_name'):
-                    vlname = varzobj.getncattr('long_name')
-                else:
-                    vlname = variables_values(zn)[4].replace('|',' ')
-                if searchInlist(varzobj.ncattrs(),'units'):
-                    vunits = varzobj.getncattr('units')
-                else:
-                    vunits = variables_values(zn)[5].replace('|',' ')
-
-                newattr = basicvardef(newvar, vsname, vlname, vunits)
-                if len(varzobj.shape) == 1:
-                    newvar[:] = varzobj[:]
-                elif len(varzobj.shape) == 2:
-                    newvar[:] = varzobj[0,:]
-                elif len(varzobj.shape) == 3:
-                    newvar[:] = varzobj[0,0,:]
-
             varvals = np.ones(tuple([Ttraj,dimz,boxs,boxs]), dtype=np.float)
+            lonvals = np.ones(tuple([Ttraj,dimz,boxs,boxs]), dtype=np.float)
+            latvals = np.ones(tuple([Ttraj,boxs,boxs]), dtype=np.float)
             rvarvals = np.ones(tuple([Ttraj,dimz,Nrad*2+1,Nrad*2+1]), dtype=np.float)
-
-            statvarvals = np.ones(tuple([Ttraj,6]), dtype=np.float)
-            rstatvarvals = np.ones(tuple([Ttraj,6]), dtype=np.float)
-
-            for it in range(dimt):
-
+            rlonvals = np.ones(tuple([Ttraj,dimz,Nrad*2+1,Nrad*2+1]), dtype=np.float)
+            rlatvals = np.ones(tuple([Ttraj,dimz,Nrad*2+1,Nrad*2+1]), dtype=np.float)
+
+            statvarvals = np.ones(tuple([Ttraj,dimz,6]), dtype=np.float)
+            rstatvarvals = np.ones(tuple([Ttraj,dimz,6]), dtype=np.float)
+
+            for it0 in range(Ttraj):
+                it = Tbeg + it0
+
+                slicev = []
+                slice2D = []
+                slicevnoT = []
+                cslicev = []
+                cslice2D = []
+                cslicevnoT = []
+
+                slicev.append(gtrajvals[it,0])
                 if gtrajvals[it,2]-box2 < 0 or gtrajvals[it,2]+box2 + 1 > dimy + 1   \
                   or gtrajvals[it,1]-box2 < 0 or gtrajvals[it,1]+box2 + 1 > dimx + 1:
-
-                    varvalst = np.ones((dimz,Nrad*2+1,Nrad*2+1),dtype=np.float)*    \
+# box values
+                    slicev.append(slice(0,dimz)
+                    slicev.append(slice(yrangeslice[it][0],yrangeslice[it][1]))
+                    slicev.append(slice(xrangeslice[it][0],xrangeslice[it][1]))
+
+                    slicevnoT.append(slice(0,dimz)
+                    slicevnoT.append(slice(yrangeslice[it][0],yrangeslice[it][1]))
+                    slicevnoT.append(slice(xrangeslice[it][0],xrangeslice[it][1]))
+
+                    slice2D.append(slice(0,dimz)
+                    slice2D.append(slice(0,yrangeslice[it][1]-yrangeslice[it][0]))
+                    slice2D.append(slice(0,xrangeslice[it][1]-xrangeslice[it][0]))
+
+                    rvarvalst = np.ones((dimz, Nrad*2+1, Nrad*2+1),dtype=np.float)*  \
                       fillValue
-
-# box values
-                    slicev.append(slice(0,dimz))
-                    slicev.append(slice(yrangeslice[it]))
-                    slicev.append(slice(xrangeslice[it]))
-
-                    slicevnoT.append(slice(0,dimz))
-                    slicevnoT.append(slice(yrangeslice[it]))
-                    slicevnoT.append(slice(xrangeslice[it]))
-
-                    slice2D.append(slice(0,dimz))
-                    slice2D.append(slice(yrangeslice2D[it]))
-                    slice2D.append(slice(xrangeslice2D[it]))
 
                     varvalst[tuple(slice2D)] = varobj[tuple(slicev)]
@@ -10608,16 +10594,18 @@
 
                 else:
+                    slicev.append(slice(0,dimz)
                     slicev.append(slice(gtrajvals[it,2]-box2, gtrajvals[it,2]+box2+1))
                     slicev.append(slice(gtrajvals[it,1]-box2, gtrajvals[it,1]+box2+1))
+                    slicevnoT.append(slice(0,dimz)
                     slicevnoT.append(slice(gtrajvals[it,2]-box2, gtrajvals[it,2]+    \
                       box2+1))
                     slicevnoT.append(slice(gtrajvals[it,1]-box2, gtrajvals[it,1]+    \
                       box2+1))
-
+                    slice2D.append(slice(0,dimz)
                     slice2D.append(slice(gtrajvals[it,2]-box2, gtrajvals[it,2]+box2 +\
                       1))
                     slice2D.append(slice(gtrajvals[it,1]-box2, gtrajvals[it,1]+box2 +\
                       1))
-                    print 'Lluis',it,'shapes varvalst:',varvals.shape,'varobj:',slicev
+
                     varvalst = varobj[tuple(slicev)]
 # box values
@@ -10637,22 +10625,23 @@
 
 # Circle values
-                cslicev.append(slice(0,dimz))
-                cslicev.append(slice(cyrangeslice[it]))
-                cslicev.append(slice(cxrangeslice[it]))
-
-                cslicevnoT.append(slice(0,dimz))
-                cslicevnoT.append(slice(cyrangeslice[it]))
-                cslicevnoT.append(slice(cxrangeslice[it]))
-
-                cslice2D.append(slice(0,dimz))
-                cslice2D.append(slice(cyrangeslice2D[it]))
-                cslice2D.append(slice(cxrangeslice2D[it]))
-
+                cslicev.append(gtrajvals[it,0])
                 if gtrajvals[it,2]-Nrad < 0 or gtrajvals[it,2]+Nrad + 1 > dimy + 1   \
                   or gtrajvals[it,1]-Nrad < 0 or gtrajvals[it,1]+Nrad + 1 > dimx + 1:
-
-                    rvarvalst = np.ones((Nrad*2+1,Nrad*2+1),dtype=np.float)*fillValue
+                    cslicev.append(slice(0,dimz)
+                    cslicev.append(slice(cyrangeslice[it][0],cyrangeslice[it][1]))
+                    cslicev.append(slice(cxrangeslice[it][0],cxrangeslice[it][1]))
+
+                    cslicevnoT.append(slice(0,dimz)
+                    cslicevnoT.append(slice(cyrangeslice[it][0],cyrangeslice[it][1]))
+                    cslicevnoT.append(slice(cxrangeslice[it][0],cxrangeslice[it][1]))
+
+                    cslice2D.append(slice(0,dimz)
+                    cslice2D.append(slice(0,cyrangeslice[it][1]-cyrangeslice[it][0]))
+                    cslice2D.append(slice(0,cxrangeslice[it][1]-cxrangeslice[it][0]))
+
+                    rvarvalst = np.ones((dimz,Nrad*2+1,Nrad*2+1),dtype=np.float)*    \
+                      fillValue
                     rvarvalst[tuple(cslice2D)] = varobj[tuple(cslicev)]
-                    rvarvalst[tuple(cslice2D)] = np.where(circdist[tuple(cslice2D)] >\
+                    rvarvalst[tuple(cslice2D)] = np.where(circdist[tuple(cslice2D)]>\
                       np.float(Nrad), fillValue, rvarvalst[tuple(cslice2D)])
 
@@ -10671,13 +10660,17 @@
 
                 else:
-                    slicev.append(slice(gtrajvals[it,2]-Nrad,gtrajvals[it,2]+Nrad+1))
-                    slicev.append(slice(gtrajvals[it,1]-Nrad,gtrajvals[it,1]+Nrad+1))
-                    slicevnoT.append(slice(gtrajvals[it,2]-Nrad, gtrajvals[it,2]+    \
+                    cslicev.append(0,dimz)
+                    cslicev.append(slice(gtrajvals[it,2]-Nrad,gtrajvals[it,2]+Nrad+1))
+                    cslicev.append(slice(gtrajvals[it,1]-Nrad,gtrajvals[it,1]+Nrad+1))
+                    cslicevnoT.append(0,dimz)
+                    cslicevnoT.append(slice(gtrajvals[it,2]-Nrad, gtrajvals[it,2]+   \
                       Nrad+1))
-                    slicevnoT.append(slice(gtrajvals[it,1]-Nrad, gtrajvals[it,1]+    \
+                    cslicevnoT.append(slice(gtrajvals[it,1]-Nrad, gtrajvals[it,1]+   \
                       Nrad+1))
-
-                    slice2D.append(slice(gtrajvals[it,2]-Nrad,gtrajvals[it,2]+Nrad+1))
-                    slice2D.append(slice(gtrajvals[it,1]-Nrad,gtrajvals[it,1]+Nrad+1))
+                    cslice2D.append(0,dimz)
+                    cslice2D.append(slice(gtrajvals[it,2]-Nrad,gtrajvals[it,2]+Nrad+ \
+                      1))
+                    cslice2D.append(slice(gtrajvals[it,1]-Nrad,gtrajvals[it,1]+Nrad+ \
+                      1))
 
                     rvarvalst = varobj[tuple(cslicev)]
@@ -10696,6 +10689,6 @@
                     maskedvals2 = maskedvals*maskedvals
                     rstatvarvals[it,:,4] = maskedvals2.mean()
-                    rstatvarvals[it,:,5] = np.sqrt(rstatvarvals[it,4] -              \
-                      rstatvarvals[it,:,3]*rstatvarvals[it,3])
+                    rstatvarvals[it,:,5] = np.sqrt(rstatvarvals[it,:,4] -            \
+                      rstatvarvals[it,:,3]*rstatvarvals[it,:,3])
 
 #            print 'statistics:',rstatvarvals[it,:] 
