me = casatools.measures()
def populate_grid(uvw, stokesI, uvgrid, uvgrid_npt, deltau, deltav, npix):
for ii in range(len(uvw[0])):
uidx = int(np.round(uvw[0][ii]//deltau + npix//2))
vidx = int(np.round(uvw[1][ii]//deltav + npix//2))
uvgrid[uidx, vidx] += stokesI[ii]
uvgrid_npt[uidx, vidx] += 1
return uvgrid, uvgrid_npt
Calculate the STD via MAD for the input data
std Calculate the std via mad
mad = np.median(np.abs(inpdat - med))
def msuvbin_to_uvgrid(ms, npix):
umin, umax = np.min(uvw[0]), np.max(uvw[0])
vmin, vmax = np.min(uvw[1]), np.max(uvw[1])
wmin, wmax = np.min(uvw[2]), np.max(uvw[2])
uvals = np.linspace(umin, umax, npix)
vvals = np.linspace(vmin, vmax, npix)
deltau = uvals[1] - uvals[0]
deltav = vvals[1] - vvals[0]
uvgrid = np.zeros((npix, npix), dtype=np.complex128)
uvgrid_npt = np.zeros((npix, npix), dtype=int)
stokesI = 0.5*(data[0] + data[1])
stokesI = np.squeeze(stokesI)
uvgrid, uvgrid_npt = populate_grid(uvw, stokesI, uvgrid, uvgrid_npt, deltau, deltav, npix)
idx = np.where(uvgrid_npt != 0)
uvgrid[idx] = uvgrid[idx]/uvgrid_npt[idx]
return uvgrid, uvgrid_npt
def calc_radial_profile(uvgrid, deltau, deltav, doplot=False):
Calculate the annular average of the uvgrid for every radius,
and fit the 1D radial profile with a polynomial.
npixx, npixy = uvgrid.shape[0], uvgrid.shape[1]
cx, cy = npixx//2, npixy//2
print(f"npixx {npixx}, npixy {npixy}, centerx {cx} centery {cy}")
x = np.arange(npixx) - cx
y = np.arange(npixy) - cy
rad = np.sqrt(x**2 + y**2).astype(int)