Commit 7bb26999 authored by Martin Reinecke's avatar Martin Reinecke

more tests, cleanups

parent 15db8450
Pipeline #51040 passed with stages
in 17 minutes and 24 seconds
......@@ -36,23 +36,23 @@ class GridderMaker(object):
raise ValueError("second dimension of uv must have length 2")
dstx, dsty = dirty_domain[0].distances
# wasteful hack to adjust to shape required by nifty_gridder
uvw = np.empty((uv.shape[0],3), dtype=np.float64)
uvw[:,0:2] = uv
uvw[:,2] = 0.
uvw = np.empty((uv.shape[0], 3), dtype=np.float64)
uvw[:, 0:2] = uv
uvw[:, 2] = 0.
# Scale uv such that 0<uv<=1 which is assumed by nifty_gridder
uvw[:, 0] = uvw[:,0]*dstx
uvw[:, 1] = uvw[:,1]*dsty
uvw[:, 0] = uvw[:, 0]*dstx
uvw[:, 1] = uvw[:, 1]*dsty
speedOfLight = 299792458.
bl = nifty_gridder.Baselines(uvw, np.array([speedOfLight]))
nxdirty, nydirty = dirty_domain.shape
nxd, nyd = dirty_domain.shape
gconf = nifty_gridder.GridderConfig(nxdirty, nydirty, eps, 1., 1.)
nu, nv = gconf.Nu(), gconf.Nv()
self._idx = nifty_gridder.getIndices(
bl, gconf, np.zeros((uv.shape[0],1),dtype=np.bool))
bl, gconf, np.zeros((uv.shape[0], 1), dtype=np.bool))
self._bl = bl
grid_domain = RGSpace([nu, nv], distances=[1, 1], harmonic=False)
du, dv = 1./(nu*dstx), 1./(nv*dsty)
grid_domain = RGSpace([nu, nv], distances=[du, dv], harmonic=True)
self._rest = _RestOperator(dirty_domain, grid_domain, gconf)
self._gridder = RadioGridder(grid_domain, bl, gconf, self._idx)
......
......@@ -41,8 +41,8 @@ def test_gridding(nu, nv, N, eps):
# Nifty
dom = ift.RGSpace((nu, nv), distances=(0.2, 1.12))
dstx, dsty = dom.distances
uv[:,0] = uv[:,0]/dstx
uv[:,1] = uv[:,1]/dsty
uv[:, 0] = uv[:, 0]/dstx
uv[:, 1] = uv[:, 1]/dsty
GM = ift.GridderMaker(dom, uv=uv, eps=eps)
vis2 = ift.from_global_data(ift.UnstructuredDomain(vis.shape), vis)
......@@ -53,10 +53,39 @@ def test_gridding(nu, nv, N, eps):
*[-ss/2 + np.arange(ss) for ss in [nu, nv]], indexing='ij')
dft = pynu*0.
for i in range(N):
dft += (vis[i]*np.exp(2j*np.pi*(x*uv[i, 0]*dstx + y*uv[i, 1]*dsty))).real
dft += (
vis[i]*np.exp(2j*np.pi*(x*uv[i, 0]*dstx + y*uv[i, 1]*dsty))).real
assert_(_l2error(dft, pynu) < eps)
def test_cartesian():
nx, ny = 2, 6
dstx, dsty = 0.3, 0.2
dom = ift.RGSpace((nx, ny), (dstx, dsty))
kx = np.fft.fftfreq(nx, dstx)
ky = np.fft.fftfreq(ny, dsty)
uu, vv = np.meshgrid(kx, ky)
tmp = np.vstack([uu[None, :], vv[None, :]])
uv = np.transpose(tmp, (2, 1, 0)).reshape(-1, 2)
GM = ift.GridderMaker(dom, uv=uv)
op = GM.getFull().adjoint
fld = ift.from_random('normal', dom)
arr = fld.to_global_data()
fld2 = ift.from_global_data(dom, np.roll(arr, (nx//2, ny//2), axis=(0, 1)))
res = op(fld2).to_global_data().reshape(nx, ny)
fft = ift.FFTOperator(dom.get_default_codomain(), target=dom).adjoint
vol = ift.full(dom, 1.).integrate()
res1 = fft(fld).to_global_data()
# FIXME: we don't understand the conjugate() yet
np.testing.assert_allclose(res, res1.conjugate()*vol)
@pmp('eps', [1e-2, 1e-6, 2e-13])
@pmp('nu', [12, 128])
@pmp('nv', [4, 12, 128])
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment