wassname/simpeg
1
0
Fork 0
mirror of https://github.com/wassname/simpeg.git synced 2026-09-10 12:37:30 +08:00
Code Issues Packages Projects Releases Wiki Activity

Testing Functions for interpolation.

Browse Source
This commit is contained in:
Rowan Cockett
2013-11-04 16:05:16 -08:00
parent df329c4955
commit 592490844d
2 changed files with 94 additions and 13 deletions
Download Patch File Download Diff File Expand all files Collapse all files
SimPEG/utils/interputils.py
+13 -13
View File
@@ -7,22 +7,22 @@ def interpmat(x,y,z,xr,yr,zr):
""" Local interpolation computed for each receiver point in turn """
nx = max(x.shape)
ny = max(y.shape)
nz = max(z.shape)
npts = max(xr.shape)
nx = x.size
ny = y.size
nz = z.size
npts = xr.shape[0]
Q = sp.lil_matrix((npts, nx*ny*nz))
for i in range(npts):
# in x-direction
# in x-direction
im = np.argmin(abs(x-xr[i]))
if xr[i] - x[im] >= 0: # Point on the left
ind_x1 = im
ind_x2 = im+1
elif xr[i] - x[im] < 0: # Point on the right
ind_x1 = im-1
ind_x2 = im
ind_x2 = im
dx1 = xr[i] - x[ind_x1]
dx2 = x[ind_x2] - xr[i]
# in y-direction
@@ -32,9 +32,9 @@ def interpmat(x,y,z,xr,yr,zr):
ind_y2 = im+1
elif yr[i] - y[im] < 0: # Point on the right
ind_y1 = im-1
ind_y2 = im
ind_y2 = im
dy1 = yr[i] - y[ind_y1]
dy2 = y[ind_y2] - yr[i]
dy2 = y[ind_y2] - yr[i]
# in z-direction
im = np.argmin(abs(z-zr[i]))
if zr[i] - z[im] >= 0: # Point on the left
@@ -42,9 +42,9 @@ def interpmat(x,y,z,xr,yr,zr):
ind_z2 = im+1
elif zr[i] - z[im] < 0: # Point on the right
ind_z1 = im-1
ind_z2 = im
ind_z2 = im
dz1 = zr[i] - z[ind_z1]
dz2 = z[ind_z2] - zr[i]
dz2 = z[ind_z2] - zr[i]
dv = (x[ind_x2] - x[ind_x1]) * (y[ind_y2] - y[ind_y1]) *(z[ind_z2] - z[ind_z1])
Dx = x[ind_x2] - x[ind_x1]
@@ -53,7 +53,7 @@ def interpmat(x,y,z,xr,yr,zr):
# Get the row in the matrix
inds = sub2ind((nx,ny,nz),[
inds = sub2ind((nx,ny,nz),[
( ind_x1, ind_y2, ind_z1),
( ind_x1, ind_y1, ind_z1),
( ind_x2, ind_y1, ind_z1),
@@ -62,7 +62,7 @@ def interpmat(x,y,z,xr,yr,zr):
( ind_x1, ind_y2, ind_z2),
( ind_x2, ind_y1, ind_z2),
( ind_x2, ind_y2, ind_z2)])
vals = [(1-dx1/Dx)*(1-dy2/Dy)*(1-dz1/Dz),
(1-dx1/Dx)*(1-dy1/Dy)*(1-dz1/Dz),
(1-dx2/Dx)*(1-dy1/Dy)*(1-dz1/Dz),
@@ -74,4 +74,4 @@ def interpmat(x,y,z,xr,yr,zr):
Q[i, mkvc(inds)] = vals
Q = Q.tocsr()
return Q
return Q