Commit fb5509a1 authored by Martin Reinecke's avatar Martin Reinecke
Browse files

remove some more over-engineering

parent 03de01dd
......@@ -50,36 +50,6 @@ void myassert(bool cond, const char *msg)
throw runtime_error(msg);
}
/*! Returns the largest integer \a n that fulfills \a 2^n<=arg. */
template<typename I> inline int ilog2 (I arg)
{
#ifdef __GNUC__
if (arg==0) return 0;
if (sizeof(I)==sizeof(int))
return 8*sizeof(int)-1-__builtin_clz(arg);
if (sizeof(I)==sizeof(long))
return 8*sizeof(long)-1-__builtin_clzl(arg);
if (sizeof(I)==sizeof(long long))
return 8*sizeof(long long)-1-__builtin_clzll(arg);
#endif
int res=0;
while (arg > 0xFFFF) { res+=16; arg>>=16; }
if (arg > 0x00FF) { res|=8; arg>>=8; }
if (arg > 0x000F) { res|=4; arg>>=4; }
if (arg > 0x0003) { res|=2; arg>>=2; }
if (arg > 0x0001) { res|=1; }
return res;
}
/*! Returns the number of bits needed to represent \a arg different values.
\a arg must be >=1. */
template<typename I> inline int bits_needed (I arg)
{
myassert(arg>=1, "argument must be >=1");
if (arg==1) return 0;
return ilog2(arg-1)+1;
}
/*! Returns the remainder of the division \a v1/v2.
The result is non-negative.
\a v1 can be positive or negative; \a v2 must be positive. */
......@@ -276,15 +246,6 @@ vector<double> correction_factors (size_t n, size_t nval, size_t w)
return res;
}
template<typename T> struct UV
{
T u, v;
UV () {}
UV (T u_, T v_) : u(u_), v(v_) {}
UV operator* (T fct) const
{ return UV(u*fct, v*fct); }
};
template<typename T> struct UVW
{
T u, v, w;
......@@ -300,7 +261,6 @@ template<typename T> class Baselines
vector<UVW<T>> coord;
vector<T> scaling;
size_t nrows, nchan;
size_t channelbits, channelmask;
public:
Baselines(const pyarr_c<T> &coord_, const pyarr_c<T> &scaling_)
......@@ -317,24 +277,19 @@ template<typename T> class Baselines
auto cood = coord_.data();
for (size_t i=0; i<coord.size(); ++i)
coord[i] = UVW<T>(cood[3*i], cood[3*i+1], cood[3*i+2]);
channelbits = bits_needed(nchan);
channelmask = (size_t(1)<<channelbits)-1;
auto rowbits = bits_needed(nrows);
myassert(rowbits+channelbits<=8*sizeof(uint32_t), "Ti too small");
myassert(nrows*nchan<(size_t(1)<<32), "too many entries in MS");
}
UVW<T> effectiveCoord(uint32_t index) const
{ return coord[index>>channelbits]*scaling[index&channelmask]; }
{
size_t irow = index/nchan;
size_t ichan = index-nchan*irow;
return coord[irow]*scaling[ichan];
}
UVW<T> effectiveCoord(size_t irow, size_t ichan) const
{ return coord[irow]*scaling[ichan]; }
size_t Nrows() const { return nrows; }
size_t Nchannels() const { return nchan; }
size_t irow(uint32_t index) const
{ return index>>channelbits; }
size_t ichannel(uint32_t index) const
{ return index&channelmask; }
size_t offset(uint32_t index) const
{ return (index>>channelbits)*nchan + (index&channelmask); }
uint32_t Index(size_t irow, size_t ichan) const
{ return (irow<<channelbits) + ichan; }
template<typename T2> pyarr_c<T2> ms2vis(const pyarr_c<T2> &ms_,
const pyarr_c<uint32_t> &idx_) const
......@@ -351,7 +306,7 @@ template<typename T> class Baselines
auto vis = res.mutable_data();
#pragma omp parallel for
for (size_t i=0; i<nvis; ++i)
vis[i] = ms[offset(idx[i])];
vis[i] = ms[idx[i]];
return res;
}
......@@ -371,7 +326,7 @@ template<typename T> class Baselines
ms[i] = T2(0);
#pragma omp parallel for
for (size_t i=0; i<nvis; ++i)
ms[offset(idx[i])] = vis[i];
ms[idx[i]] = vis[i];
return res;
}
......@@ -391,7 +346,7 @@ template<typename T> class Baselines
#pragma omp parallel for
for (size_t i=0; i<nvis; ++i)
ms[offset(idx[i])] += vis[i];
ms[idx[i]] += vis[i];
return ms_;
}
};
......@@ -916,8 +871,7 @@ template<typename T> pyarr_c<uint32_t> getIndices(const Baselines<T> &baselines,
for (int ichan=chbegin; ichan<chend; ++ichan)
if (!flags[irow*nchan + ichan])
{
auto idx = baselines.Index(irow, ichan);
auto uvw = baselines.effectiveCoord(idx);
auto uvw = baselines.effectiveCoord(irow, ichan);
if ((uvw.w>=wmin) && (uvw.w<wmax))
{
T u, v;
......@@ -938,8 +892,7 @@ template<typename T> pyarr_c<uint32_t> getIndices(const Baselines<T> &baselines,
for (int ichan=chbegin; ichan<chend; ++ichan)
if (!flags[irow*nchan + ichan])
{
auto idx = baselines.Index(irow, ichan);
auto uvw = baselines.effectiveCoord(idx);
auto uvw = baselines.effectiveCoord(irow, ichan);
if ((uvw.w>=wmin) && (uvw.w<wmax))
{
T u, v;
......@@ -947,7 +900,7 @@ template<typename T> pyarr_c<uint32_t> getIndices(const Baselines<T> &baselines,
gconf.getpix(uvw.u, uvw.v, u, v, iu0, iv0);
iu0 = (iu0+nsafe)>>logsquare;
iv0 = (iv0+nsafe)>>logsquare;
iout[acc[nbv*iu0 + iv0]] = idx;
iout[acc[nbv*iu0 + iv0]] = irow*nchan + ichan;
++acc[nbv*iu0 + iv0];
}
}
......
Supports Markdown
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