Commit f3ba8494 authored by Martin Reinecke's avatar Martin Reinecke

Merge branch 'ducc0' into mpi_experiments

parents cc802abd 9ce68e92
0.7.0:
- wgridder:
- performance (especially scaling) improvements
- oversampling factors up to 2.5 supported
- new, more flexible interface in submodule `wgridder.experimental`
(subject to further changes!)
0.6.0:
- general:
- multi-threading improvements contributed by Peter Bell
......@@ -21,7 +29,6 @@
0.2.0:
- wgridder:
- kernels are now evaluated via polynomial approximation, allowing much
more freedom in the choice of kernel function
......
[build-system]
requires = ["setuptools >= 40.6.0", "pybind11 >= 2.5.0", "numpy >= 1.17.0"]
requires = ["setuptools >= 40.6.0", "pybind11 >= 2.6.0", "numpy >= 1.17.0"]
build-backend = "setuptools.build_meta"
......@@ -15,7 +15,7 @@
from time import time
import ducc0.wgridder as wgridder
import ducc0.wgridder.experimental as wgridder
import matplotlib.pyplot as plt
import numpy as np
......@@ -36,21 +36,22 @@ def main():
pixsize = fov_deg/npixdirty*DEG2RAD
nthreads = 2
epsilon = 1e-4
print('Start gridding...')
do_wstacking = True
do_wgridding = True
print('Start gridding...')
t0 = time()
dirty = wgridder.ms2dirty(uvw, freq, vis, wgt, npixdirty, npixdirty, pixsize,
pixsize, 0, 0, epsilon, do_wstacking, nthreads, verbosity=1, mask=flags)
print('Done')
dirty = wgridder.vis2dirty(uvw=uvw, freq=freq, vis=vis, wgt=wgt,
npix_x=npixdirty, npix_y=npixdirty, pixsize_x=pixsize,
pixsize_y=pixsize, epsilon=epsilon, do_wgridding=do_wgridding,
nthreads=nthreads, verbosity=1, mask=flags, flip_v=True)
t = time() - t0
print("{} s".format(t))
print("{} visibilities/thread/s".format(np.sum(wgt != 0)/nthreads/t))
t0 = time()
_ = wgridder.dirty2ms(uvw, freq, dirty, wgt, pixsize,
pixsize, 0, 0, epsilon, do_wstacking, nthreads, verbosity=1, mask=flags)
print('Done')
wgridder.dirty2vis(uvw=uvw, freq=freq, dirty=dirty, wgt=wgt,
pixsize_x=pixsize, pixsize_y=pixsize, epsilon=epsilon,
do_wgridding=do_wgridding, nthreads=nthreads, verbosity=1, mask=flags,
flip_v=True)
t = time() - t0
print("{} s".format(t))
print("{} visibilities/thread/s".format(np.sum(wgt != 0)/nthreads/t))
......
......@@ -673,7 +673,7 @@ out : int
} // unnamed namespace
void add_fft(py::module &msup)
void add_fft(py::module_ &msup)
{
using namespace pybind11::literals;
auto m = msup.def_submodule("fft");
......
......@@ -73,7 +73,7 @@ template<typename T> void quickzero(mav<T,2> &arr, size_t nthreads)
}
template<typename T> complex<T> hsum_cmplx(native_simd<T> vr, native_simd<T> vi)
{ return complex<T>(reduce(vr, std::plus<>()), reduce(vi, std::plus<>())); }
{ return complex<T>(reduce(vr, plus<>()), reduce(vi, plus<>())); }
#if (defined(__AVX__) && (!defined(__AVX512F__)))
#if 1
......@@ -147,7 +147,7 @@ template<typename T> void hartley2complex
size_t xv = (v==0) ? 0 : nv-v;
T v1 = T(0.5)*grid( u, v);
T v2 = T(0.5)*grid(xu,xv);
grid2.v(u,v) = std::complex<T>(v1+v2, v1-v2);
grid2.v(u,v) = complex<T>(v1+v2, v1-v2);
}
}
});
......@@ -274,6 +274,8 @@ class Baselines
UVW effectiveCoord(size_t row, size_t chan) const
{ return coord[row]*f_over_c[chan]; }
double absEffectiveW(size_t row, size_t chan) const
{ return abs(coord[row].w*f_over_c[chan]); }
UVW baseCoord(size_t row) const
{ return coord[row]; }
double ffact(size_t chan) const
......@@ -305,6 +307,7 @@ template<typename T> class Params
size_t nthreads;
size_t verbosity;
bool negate_v, divide_by_n;
double sigma_min, sigma_max;
Baselines bl;
VVR ranges;
......@@ -570,7 +573,7 @@ template<typename T> class Params
void countRanges()
{
timers.push("range count");
timers.push("building index");
size_t nrow=bl.Nrows(),
nchan=bl.Nchannels();
......@@ -592,142 +595,108 @@ template<typename T> class Params
constexpr double max_asymm = 0.01;
size_t max_allowed = size_t(nvis/double(nbunch*nthreads)*max_asymm);
struct tmp1
{
vector<RowchanRange> v;
size_t sz=0;
tmp1()
: sz(0)
{ v.reserve(8); }
void add(const RowchanRange &rng)
{
v.push_back(rng);
sz += rng.ch_end-rng.ch_begin;
}
};
struct tmp2
{
vector<tmp1> v;
size_t sz=0;
vector<vector<RowchanRange>> v;
void add(const RowchanRange &rng, size_t max_allowed)
{
if (v.empty() || (v.back().sz>=max_allowed))
v.push_back(tmp1());
v.back().add(rng);
}
void add(tmp2 &&other, size_t max_allowed)
{
if (other.v.empty()) return;
if (v.empty())
{ v=other.v; return; }
auto backup = v.back();
v.pop_back();
v.reserve(v.size()+other.v.size()+2);
for (auto &&x: other.v)
v.push_back(x);
for (auto &&x: backup.v)
add(x, max_allowed);
if (v.empty() || (sz>=max_allowed))
{
v.emplace_back();
sz=0;
}
v.back().push_back(rng);
sz += rng.ch_end-rng.ch_begin;
}
};
using Vmap = map<Uvwidx, tmp2>;
struct bufmap
{
Vmap m;
mutex mut;
uint64_t dummy[8]; // separator to keep every entry on a different cache line
};
vector<bufmap> buf(nthreads);
bool have_wgt=wgt.size()!=0;
if (have_wgt) checkShape(wgt.shape(),{nrow,nchan});
bool have_ms=ms_in.size()!=0;
if (have_ms) checkShape(ms_in.shape(), {nrow,nchan});
bool have_mask=mask.size()!=0;
if (have_mask) checkShape(mask.shape(), {nrow,nchan});
checkShape(wgt.shape(),{nrow,nchan});
checkShape(ms_in.shape(), {nrow,nchan});
checkShape(mask.shape(), {nrow,nchan});
size_t ntiles_u = (nu>>logsquare) + 20,
ntiles_v = (nv>>logsquare) + 20;
vector<bufmap> buf(ntiles_u*ntiles_v);
auto chunk = max<size_t>(1, nrow/(20*nthreads));
auto xdw = 1./dw;
auto shift = dw-(0.5*supp*dw)-wmin;
execDynamic(nrow, nthreads, chunk, [&](Scheduler &sched)
{
auto &mymap(buf[sched.thread_num()].m);
vector<pair<uint16_t, uint16_t>> interbuf;
while (auto rng=sched.getNext())
for(auto irow=rng.lo; irow<rng.hi; ++irow)
{
bool on=false;
Uvwidx uvwlast(0,0,0);
tmp2 *ptr=0;
size_t chan0=0;
auto flush=[&]()
{
if (interbuf.empty()) return;
auto tileidx = uvwlast.tile_u + ntiles_u*uvwlast.tile_v;
lock_guard<mutex> lock(buf[tileidx].mut);
auto &loc(buf[tileidx].m[uvwlast]);
for (auto &x: interbuf)
loc.add(RowchanRange(irow, x.first, x.second), max_allowed);
interbuf.clear();
};
auto add=[&](uint16_t cb, uint16_t ce)
{ interbuf.emplace_back(cb, ce); };
for (size_t ichan=0; ichan<nchan; ++ichan)
{
if (((!have_ms ) || (norm(ms_in(irow,ichan))!=0)) &&
((!have_wgt) || (wgt(irow,ichan)!=0)) &&
((!have_mask) || (mask(irow,ichan)!=0)))
if (norm(ms_in(irow,ichan))*wgt(irow,ichan)*mask(irow,ichan)!=0)
{
auto uvw = bl.effectiveCoord(irow, ichan);
if (uvw.w<0) uvw.Flip();
uvw.FixW();
double u, v;
int iu0, iv0, iw;
getpix(uvw.u, uvw.v, u, v, iu0, iv0);
iu0 = (iu0+nsafe)>>logsquare;
iv0 = (iv0+nsafe)>>logsquare;
iw = do_wgridding ?
max(0,int(1+(abs(uvw.w)-(0.5*supp*dw)-wmin)/dw)) : 0;
iw = do_wgridding ? max(0,int((uvw.w+shift)*xdw)) : 0;
Uvwidx uvwcur(iu0, iv0, iw);
if (!on) // new active region
{
on=true;
if ((!ptr) || (uvwlast!=uvwcur)) ptr=&mymap[uvwcur];
uvwlast = uvwcur;
chan0=ichan;
if (uvwlast!=uvwcur) flush();
uvwlast=uvwcur; chan0=ichan;
}
else if (uvwlast!=uvwcur) // change of active region
{
ptr->add(RowchanRange(irow, chan0, ichan), max_allowed);
add(chan0, ichan);
flush();
uvwlast=uvwcur; chan0=ichan;
ptr=&mymap[uvwcur];
}
}
else if (on) // end of active region
{
ptr->add(RowchanRange(irow, chan0, ichan), max_allowed);
add(chan0, ichan);
on=false;
}
}
if (on) // end of active region at last channel
ptr->add(RowchanRange(irow, chan0, nchan), max_allowed);
add(chan0, nchan);
flush();
}
});
timers.poppush("range merging");
size_t nth = nthreads;
while (nth>1)
{
auto nmerge=nth/2;
execParallel(nmerge, [&](Scheduler &sched)
{
auto tid = sched.thread_num();
auto &s1 = buf[tid].m;
auto &s2 = buf[nth-1-tid].m;
for (auto &&v : s2)
{
auto loc = s1.find(v.first);
if (loc == s1.end())
s1[v.first] = move(v.second);
else
loc->second.add(move(v.second), max_allowed);
}
Vmap().swap(s2);
});
nth-=nmerge;
}
timers.poppush("building final range vector");
size_t total=0;
for (const auto &x: buf[0].m)
total += x.second.v.size();
for (const auto &x: buf)
for (const auto &y: x.m)
total += y.second.v.size();
ranges.reserve(total);
for (auto &v : buf[0].m)
{
for (auto &v2:v.second.v)
ranges.emplace_back(v.first, move(v2.v));
}
for (auto &x: buf)
for (auto &y: x.m)
for (auto &z: y.second.v)
ranges.emplace_back(y.first, move(z));
timers.pop();
}
......@@ -752,7 +721,7 @@ template<typename T> class Params
mav<T,2> bufr, bufi;
T *px0r, *px0i;
double w0, xdw;
vector<std::mutex> &locks;
vector<mutex> &locks;
DUCC0_NOINLINE void dump()
{
......@@ -766,7 +735,7 @@ template<typename T> class Params
{
int idxv = idxv0;
{
std::lock_guard<std::mutex> lock(locks[idxu]);
lock_guard<mutex> lock(locks[idxu]);
for (int iv=0; iv<sv; ++iv)
{
grid.v(idxu,idxv) += complex<T>(bufr(iu,iv), bufi(iu,iv));
......@@ -787,7 +756,7 @@ template<typename T> class Params
kbuf buf;
HelperX2g2(const Params *parent_, mav<complex<T>,2> &grid_,
vector<std::mutex> &locks_, double w0_=-1, double dw_=-1)
vector<mutex> &locks_, double w0_=-1, double dw_=-1)
: parent(parent_), tkrn(*parent->krn), grid(grid_),
iu0(-1000000), iv0(-1000000),
bu0(-1000000), bv0(-1000000),
......@@ -920,8 +889,7 @@ template<typename T> class Params
template<size_t SUPP, bool wgrid> [[gnu::hot]] void x2grid_c_helper
(mav<complex<T>,2> &grid, size_t p0, double w0)
{
bool have_wgt = wgt.size()!=0;
vector<std::mutex> locks(nu);
vector<mutex> locks(nu);
execDynamic(ranges.size(), nthreads, wgrid ? SUPP : 1, [&](Scheduler &sched)
{
......@@ -950,7 +918,7 @@ auto ix = ix_+ranges.size()/2; if (ix>=ranges.size()) ix -=ranges.size();
hlp.prep(coord, nth);
auto v(ms_in(row, ch));
if (have_wgt) v*=wgt(row, ch);
v*=wgt(row, ch);
if constexpr (NVEC==1)
{
......@@ -1033,8 +1001,6 @@ auto ix = ix_+ranges.size()/2; if (ix>=ranges.size()) ix -=ranges.size();
template<size_t SUPP, bool wgrid> [[gnu::hot]] void grid2x_c_helper
(const mav<complex<T>,2> &grid, size_t p0, double w0)
{
bool have_wgt = wgt.size()!=0;
// Loop over sampling points
execDynamic(ranges.size(), nthreads, wgrid ? SUPP : 1, [&](Scheduler &sched)
{
......@@ -1092,7 +1058,7 @@ auto ix = ix_+ranges.size()/2; if (ix>=ranges.size()) ix -=ranges.size();
}
ri *= imflip;
auto r = hsum_cmplx(rr,ri);
if (have_wgt) r*=wgt(row, ch);
r*=wgt(row, ch);
ms_out.v(row, ch) += r;
}
}
......@@ -1300,7 +1266,7 @@ auto ix = ix_+ranges.size()/2; if (ix>=ranges.size()) ix -=ranges.size();
nm1min = sqrt(max(1.-x0*x0-y0*y0,0.))-1.;
if (x0*x0+y0*y0>1.)
nm1min = -sqrt(abs(1.-x0*x0-y0*y0))-1.;
auto idx = getAvailableKernels<T>(epsilon);
auto idx = getAvailableKernels<T>(epsilon, sigma_min, sigma_max);
double mincost = 1e300;
constexpr double nref_fft=2048;
constexpr double costref_fft=0.0693;
......@@ -1344,12 +1310,9 @@ auto ix = ix_+ranges.size()/2; if (ix>=ranges.size()) ix -=ranges.size();
timers.push("Initial scan");
size_t nrow=bl.Nrows(),
nchan=bl.Nchannels();
bool have_wgt=wgt.size()!=0;
if (have_wgt) checkShape(wgt.shape(),{nrow,nchan});
bool have_ms=ms_in.size()!=0;
if (have_ms) checkShape(ms_in.shape(), {nrow,nchan});
bool have_mask=mask.size()!=0;
if (have_mask) checkShape(mask.shape(), {nrow,nchan});
checkShape(wgt.shape(),{nrow,nchan});
checkShape(ms_in.shape(), {nrow,nchan});
checkShape(mask.shape(), {nrow,nchan});
nvis=0;
wmin_d=1e300;
......@@ -1361,13 +1324,10 @@ auto ix = ix_+ranges.size()/2; if (ix>=ranges.size()) ix -=ranges.size();
size_t lnvis=0;
for(auto irow=lo; irow<hi; ++irow)
for (size_t ichan=0, idx=irow*nchan; ichan<nchan; ++ichan, ++idx)
if (((!have_ms ) || (norm(ms_in(irow,ichan))!=0)) &&
((!have_wgt) || (wgt(irow,ichan)!=0)) &&
((!have_mask) || (mask(irow,ichan)!=0)))
if (norm(ms_in(irow,ichan))*wgt(irow,ichan)*mask(irow,ichan) != 0)
{
++lnvis;
auto uvw = bl.effectiveCoord(irow,ichan);
double w = abs(uvw.w);
double w = bl.absEffectiveW(irow, ichan);
lwmin_d = min(lwmin_d, w);
lwmax_d = max(lwmax_d, w);
}
......@@ -1388,8 +1348,9 @@ auto ix = ix_+ranges.size()/2; if (ix>=ranges.size()) ix -=ranges.size();
const mav<T,2> &wgt_, const mav<uint8_t,2> &mask_,
double pixsize_x_, double pixsize_y_, double epsilon_,
bool do_wgridding_, size_t nthreads_, size_t verbosity_,
bool negate_v_, bool divide_by_n_)
: gridding(ms_in_.size()>0),
bool negate_v_, bool divide_by_n_, double sigma_min_,
double sigma_max_)
: gridding(ms_out_.size()==0),
timers(gridding ? "gridding" : "degridding"),
ms_in(ms_in_), ms_out(ms_out_),
dirty_in(dirty_in_), dirty_out(dirty_out_),
......@@ -1400,7 +1361,8 @@ auto ix = ix_+ranges.size()/2; if (ix>=ranges.size()) ix -=ranges.size();
epsilon(epsilon_),
do_wgridding(do_wgridding_),
nthreads(nthreads_), verbosity(verbosity_),
negate_v(negate_v_), divide_by_n(divide_by_n_)
negate_v(negate_v_), divide_by_n(divide_by_n_),
sigma_min(sigma_min_), sigma_max(sigma_max_)
{
timers.push("Baseline construction");
bl = Baselines(uvw, freq, negate_v);
......@@ -1458,28 +1420,36 @@ auto ix = ix_+ranges.size()/2; if (ix>=ranges.size()) ix -=ranges.size();
}
};
// Note to self: divide_by_n should always be true when doing Bayesian imaging,
// but wsclean needs it to be false, so this must be kept as a parameter.
template<typename T> void ms2dirty(const mav<double,2> &uvw,
const mav<double,1> &freq, const mav<complex<T>,2> &ms,
const mav<T,2> &wgt, const mav<uint8_t,2> &mask, double pixsize_x, double pixsize_y, double epsilon,
const mav<T,2> &wgt_, const mav<uint8_t,2> &mask_, double pixsize_x, double pixsize_y, double epsilon,
bool do_wgridding, size_t nthreads, mav<T,2> &dirty, size_t verbosity,
bool negate_v=false, bool divide_by_n=true)
bool negate_v=false, bool divide_by_n=true, double sigma_min=1.1,
double sigma_max=2.6)
{
mav<complex<T>,2> ms_out(nullptr, {0,0}, false);
mav<T,2> dirty_in(nullptr, {0,0}, false);
Params<T> par(uvw, freq, ms, ms_out, dirty_in, dirty, wgt, mask, pixsize_x, pixsize_y, epsilon, do_wgridding, nthreads, verbosity, negate_v, divide_by_n);
auto ms_out(mav<complex<T>,2>::build_empty());
auto dirty_in(mav<T,2>::build_empty());
auto wgt(wgt_.size()!=0 ? wgt_ : wgt_.build_uniform(ms.shape(), 1.));
auto mask(mask_.size()!=0 ? mask_ : mask_.build_uniform(ms.shape(), 1));
Params<T> par(uvw, freq, ms, ms_out, dirty_in, dirty, wgt, mask, pixsize_x,
pixsize_y, epsilon, do_wgridding, nthreads, verbosity, negate_v,
divide_by_n, sigma_min, sigma_max);
}
template<typename T> void dirty2ms(const mav<double,2> &uvw,
const mav<double,1> &freq, const mav<T,2> &dirty,
const mav<T,2> &wgt, const mav<uint8_t,2> &mask, double pixsize_x, double pixsize_y,
const mav<T,2> &wgt_, const mav<uint8_t,2> &mask_, double pixsize_x, double pixsize_y,
double epsilon, bool do_wgridding, size_t nthreads, mav<complex<T>,2> &ms,
size_t verbosity, bool negate_v=false, bool divide_by_n=true)
size_t verbosity, bool negate_v=false, bool divide_by_n=true,
double sigma_min=1.1, double sigma_max=2.6)
{
mav<complex<T>,2> ms_in(nullptr, {0,0}, false);
mav<T,2> dirty_out(nullptr, {0,0}, false);
Params<T> par(uvw, freq, ms_in, ms, dirty, dirty_out, wgt, mask, pixsize_x, pixsize_y, epsilon, do_wgridding, nthreads, verbosity, negate_v, divide_by_n);
auto ms_in(mav<complex<T>,2>::build_uniform(ms.shape(),1.));
auto dirty_out(mav<T,2>::build_empty());
auto wgt(wgt_.size()!=0 ? wgt_ : wgt_.build_uniform(ms.shape(), 1.));
auto mask(mask_.size()!=0 ? mask_ : mask_.build_uniform(ms.shape(), 1));
Params<T> par(uvw, freq, ms_in, ms, dirty, dirty_out, wgt, mask, pixsize_x,
pixsize_y, epsilon, do_wgridding, nthreads, verbosity, negate_v,
divide_by_n, sigma_min, sigma_max);
}
} // namespace detail_gridder
......
......@@ -346,7 +346,7 @@ that their last dimension is removed.
The employed algorithm is highly accurate, even for angles close to 0 or pi.
)""";
void add_healpix(py::module &msup)
void add_healpix(py::module_ &msup)
{
using namespace pybind11::literals;
auto m = msup.def_submodule("healpix");
......
......@@ -208,7 +208,7 @@ const char *misc_DS = R"""(
Various unsorted utilities
)""";
void add_misc(py::module &msup)
void add_misc(py::module_ &msup)
{
using namespace pybind11::literals;
auto m = msup.def_submodule("misc");
......
......@@ -203,7 +203,7 @@ numpy.ndarray((nval, 4), dtype=numpy.float64) : the output quaternions
This is identical to the provided "out" array.
)""";
void add_pointingprovider(py::module &msup)
void add_pointingprovider(py::module_ &msup)
{
using namespace pybind11::literals;
auto m = msup.def_submodule("pointingprovider");
......
......@@ -189,7 +189,7 @@ Python interface for some of the libsharp functionality
Error conditions are reported by raising exceptions.
)""";
void add_sht(py::module &msup)
void add_sht(py::module_ &msup)
{
using namespace pybind11::literals;
auto m = msup.def_submodule("sht");
......
......@@ -235,7 +235,7 @@ Notes
- must be the last call to the object
)""";
void add_totalconvolve(py::module &msup)
void add_totalconvolve(py::module_ &msup)
{
using namespace pybind11::literals;
auto m = msup.def_submodule("totalconvolve");
......
......@@ -34,42 +34,209 @@ namespace py = pybind11;
auto None = py::none();
template<typename T> py::array ms2dirty2(const py::array &uvw_,
const py::array &freq_, const py::array &ms_, const py::object &wgt_, const py::object &mask_,
template<typename T> py::array vis2dirty2(const py::array &uvw_,
const py::array &freq_, const py::array &vis_, const py::object &wgt_, const py::object &mask_,
size_t npix_x, size_t npix_y, double pixsize_x, double pixsize_y,
double epsilon, bool do_wgridding, size_t nthreads,
size_t verbosity)
double epsilon, bool do_wgridding, size_t nthreads, size_t verbosity,
bool flip_v, bool divide_by_n, py::object &dirty_, double sigma_min,
double sigma_max)
{
auto uvw = to_mav<double,2>(uvw_, false);
auto freq = to_mav<double,1>(freq_, false);
auto ms = to_mav<complex<T>,2>(ms_, false);
auto wgt = get_optional_const_Pyarr<T>(wgt_, {ms.shape(0),ms.shape(1)});
auto vis = to_mav<complex<T>,2>(vis_, false);
auto wgt = get_optional_const_Pyarr<T>(wgt_, {vis.shape(0),vis.shape(1)});
auto wgt2 = to_mav<T,2>(wgt, false);
auto mask = get_optional_const_Pyarr<uint8_t>(mask_, {uvw.shape(0),freq.shape(0)});
auto mask2 = to_mav<uint8_t,2>(mask, false);
auto dirty = make_Pyarr<T>({npix_x,npix_y});
// sizes must be either both zero or both nonzero
MR_assert((npix_x==0)==(npix_y==0), "inconsistent dirty image dimensions");
auto dirty = (npix_x==0) ? get_Pyarr<T>(dirty_, 2)
: get_optional_Pyarr<T>(dirty_, {npix_x, npix_y});
auto dirty2 = to_mav<T,2>(dirty, true);
{
py::gil_scoped_release release;
ms2dirty(uvw,freq,ms,wgt2,mask2,pixsize_x,pixsize_y,epsilon,
do_wgridding,nthreads,dirty2,verbosity);
ms2dirty(uvw,freq,vis,wgt2,mask2,pixsize_x,pixsize_y,epsilon,
do_wgridding,nthreads,dirty2,verbosity,flip_v,divide_by_n, sigma_min,
sigma_max);
}
return move(dirty);
}
py::array Pyvis2dirty(const py::array &uvw,
const py::array &freq, const py::array &vis, const py::object &wgt,
size_t npix_x, size_t npix_y, double pixsize_x, double pixsize_y,
double epsilon, bool do_wgridding, size_t nthreads,
size_t verbosity, const py::object &mask, bool flip_v, bool divide_by_n,
py::object &dirty=None, double sigma_min=1.1, double sigma_max=2.6)
{
if (isPyarr<complex<float>>(vis))
return vis2dirty2<float>(uvw, freq, vis, wgt, mask, npix_x, npix_y,
pixsize_x, pixsize_y, epsilon, do_wgridding, nthreads, verbosity,
flip_v, divide_by_n, dirty, sigma_min, sigma_max);
if (isPyarr<complex<double>>(vis))
return vis2dirty2<double>(uvw, freq, vis, wgt, mask, npix_x, npix_y,
<