Commit 632fa522 authored by Martin Reinecke's avatar Martin Reinecke

fix and cleanup

parent 04c1d771
......@@ -1317,9 +1317,8 @@ template<typename T> void update_idx(vector<T> &v, vector<T> &vold,
}
template<typename Serv> void wstack_common(
const GridderConfig &gconf, const Serv &srv, double &wmin,
double &dw, size_t &nplanes, vector<size_t> &nvis_plane,
vector<vector<idx_t>> &minplane, size_t verbosity)
const GridderConfig &gconf, const Serv &srv, double &wmin, double &dw,
size_t &nplanes, vector<vector<idx_t>> &minplane, size_t verbosity)
{
size_t nvis = srv.Nvis();
size_t nthreads = gconf.Nthreads();
......@@ -1345,7 +1344,6 @@ template<typename Serv> void wstack_common(
nplanes += supp-2;
if (verbosity>0) cout << "Kernel support: " << supp << endl;
if (verbosity>0) cout << "nplanes: " << nplanes << endl;
nvis_plane.resize(nplanes,0);
minplane.resize(nplanes);
vector<size_t> tcnt(my_max_threads()*nplanes,0);
......@@ -1363,12 +1361,8 @@ template<typename Serv> void wstack_common(
++nvp[i];
}
#pragma omp critical (wstack_common)
{
for (size_t i=0; i<nplanes; ++i)
nvis_plane[i] += nvp[i];
for (size_t i=0; i<nplanes; ++i)
tcnt[mythread*nplanes+i] = mytcnt[i];
}
#pragma omp barrier
#pragma omp single
for (size_t j=0; j<nplanes; ++j)
......@@ -1405,10 +1399,9 @@ template<typename T, typename Serv> void x2dirty(
double wmin;
double dw;
size_t nplanes;
vector<size_t> nvis_plane;
vector<vector<idx_t>> minplane;
if (verbosity>0) cout << "Gridding using improved w-stacking" << endl;
wstack_common(gconf, srv, wmin, dw, nplanes, nvis_plane, minplane, verbosity);
wstack_common(gconf, srv, wmin, dw, nplanes, minplane, verbosity);
dirty.fill(0);
vector<idx_t> subidx, subidx_old;
tmpStorage<complex<T>,2> grid_({gconf.Nu(),gconf.Nv()});
......@@ -1417,13 +1410,13 @@ template<typename T, typename Serv> void x2dirty(
auto tdirty=tdirty_.getMav();
for (size_t iw=0; iw<nplanes; ++iw)
{
if (nvis_plane[iw]==0) continue;
update_idx(subidx, subidx_old, minplane[iw], iw>=supp ? minplane[iw-supp] : vector<idx_t>());
if (subidx.size()==0) continue;
if (verbosity>1)
cout << "Working on plane " << iw << " containing " << nvis_plane[iw]
cout << "Working on plane " << iw << " containing " << subidx.size()
<< " visibilities" << endl;
double wcur = wmin+iw*dw;
update_idx(subidx, subidx_old, minplane[iw], iw>=supp ? minplane[iw-supp] : vector<idx_t>());
auto subidx2 = const_mav<idx_t, 1>(subidx.data(), {subidx.size()});
auto subsrv = srv.getSubserv(subidx2);
grid.fill(0);
......@@ -1477,10 +1470,9 @@ template<typename T, typename Serv> void dirty2x(
double wmin;
double dw;
size_t nplanes;
vector<size_t> nvis_plane;
vector<vector<idx_t>> minplane;
if (verbosity>0) cout << "Degridding using improved w-stacking" << endl;
wstack_common(gconf, srv, wmin, dw, nplanes, nvis_plane, minplane, verbosity);
wstack_common(gconf, srv, wmin, dw, nplanes, minplane, verbosity);
tmpStorage<T,2> tdirty_({nx_dirty,ny_dirty});
auto tdirty=tdirty_.getMav();
......@@ -1497,9 +1489,10 @@ template<typename T, typename Serv> void dirty2x(
auto tdirty2=tdirty2_.getMav();
for (size_t iw=0; iw<nplanes; ++iw)
{
if (nvis_plane[iw]==0) continue;
update_idx(subidx, subidx_old, minplane[iw], iw>=supp ? minplane[iw-supp] : vector<idx_t>());
if (subidx.size()==0) continue;
if (verbosity>1)
cout << "Working on plane " << iw << " containing " << nvis_plane[iw]
cout << "Working on plane " << iw << " containing " << subidx.size()
<< " visibilities" << endl;
double wcur = wmin+iw*dw;
......@@ -1510,7 +1503,6 @@ template<typename T, typename Serv> void dirty2x(
gconf.apply_wscreen(cmav(tdirty2), tdirty2, wcur, false);
gconf.dirty2grid_c(cmav(tdirty2), grid);
update_idx(subidx, subidx_old, minplane[iw], iw>=supp ? minplane[iw-supp] : vector<idx_t>());
auto subidx2 = const_mav<idx_t, 1>(subidx.data(), {subidx.size()});
auto subsrv = srv.getSubserv(subidx2);
grid2x_c(gconf, cmav(grid), subsrv, wcur, dw);
......
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