PURIFY
Next-generation radio interferometric imaging
Functions
purify::measurementoperator Namespace Reference

Functions

template<class T >
std::shared_ptr< sopt::LinearTransform< T > > init_degrid_operator_2d (const Vector< t_real > &u, const Vector< t_real > &v, const Vector< t_real > &w, const Vector< t_complex > &weights, const t_uint &imsizey, const t_uint &imsizex, const t_real &oversample_ratio=2, const kernels::kernel kernel=kernels::kernel::kb, const t_uint Ju=4, const t_uint Jv=4, const bool w_stacking=false, const t_real &cellx=1, const t_real &celly=1)
 Returns linear transform that is the standard degridding operator. More...
 
template<class T >
std::shared_ptr< sopt::LinearTransform< T > > init_degrid_operator_2d (const utilities::vis_params &uv_vis_input, const t_uint &imsizey, const t_uint &imsizex, const t_real &cell_x, const t_real &cell_y, const t_real &oversample_ratio=2, const kernels::kernel kernel=kernels::kernel::kb, const t_uint Ju=4, const t_uint Jv=4, const bool w_stacking=false)
 
template<class T >
std::shared_ptr< sopt::LinearTransform< T > > init_degrid_operator_2d (const Vector< t_real > &u, const Vector< t_real > &v, const Vector< t_real > &w, const Vector< t_complex > &weights, const t_uint imsizey, const t_uint imsizex, const t_real oversample_ratio, const kernels::kernel kernel, const t_uint Ju, const t_uint Jw, const bool w_stacking, const t_real cellx, const t_real celly, const t_real absolute_error, const t_real relative_error, const dde_type dde)
 Returns linear transform that is the standard degridding operator. More...
 
template<class T >
std::shared_ptr< sopt::LinearTransform< T > > init_degrid_operator_2d (const utilities::vis_params &uv_vis_input, const t_uint imsizey, const t_uint imsizex, const t_real cell_x, const t_real cell_y, const t_real oversample_ratio, const kernels::kernel kernel, const t_uint Ju, const t_uint Jw, const bool w_stacking, const t_real absolute_error, const t_real relative_error, const dde_type dde)
 

Function Documentation

◆ init_degrid_operator_2d() [1/4]

template<class T >
std::shared_ptr<sopt::LinearTransform<T> > purify::measurementoperator::init_degrid_operator_2d ( const utilities::vis_params uv_vis_input,
const t_uint &  imsizey,
const t_uint &  imsizex,
const t_real &  cell_x,
const t_real &  cell_y,
const t_real &  oversample_ratio = 2,
const kernels::kernel  kernel = kernels::kernel::kb,
const t_uint  Ju = 4,
const t_uint  Jv = 4,
const bool  w_stacking = false 
)

Definition at line 625 of file operators.h.

629  {
630  const auto uv_vis = utilities::convert_to_pixels(uv_vis_input, cell_x, cell_y, imsizex, imsizey,
631  oversample_ratio);
632  return init_degrid_operator_2d<T>(uv_vis.u, uv_vis.v, uv_vis.w, uv_vis.weights, imsizey, imsizex,
633  oversample_ratio, kernel, Ju, Jv, w_stacking, cell_x, cell_y);
634 }
utilities::vis_params w_stacking(utilities::vis_params const &params, sopt::mpi::Communicator const &comm, const t_int iters, const std::function< t_real(t_real)> &cost, const t_real k_means_rel_diff)
utilities::vis_params convert_to_pixels(const utilities::vis_params &uv_vis, const t_real cell_x, const t_real cell_y, const t_real imsizex, const t_real imsizey, const t_real oversample_ratio)
Converts u and v coordaintes to units of pixels.

References purify::utilities::convert_to_pixels(), and purify::utilities::w_stacking().

◆ init_degrid_operator_2d() [2/4]

template<class T >
std::shared_ptr<sopt::LinearTransform<T> > purify::measurementoperator::init_degrid_operator_2d ( const utilities::vis_params uv_vis_input,
const t_uint  imsizey,
const t_uint  imsizex,
const t_real  cell_x,
const t_real  cell_y,
const t_real  oversample_ratio,
const kernels::kernel  kernel,
const t_uint  Ju,
const t_uint  Jw,
const bool  w_stacking,
const t_real  absolute_error,
const t_real  relative_error,
const dde_type  dde 
)

Definition at line 207 of file wproj_operators.h.

211  {
212  const auto uv_vis = utilities::convert_to_pixels(uv_vis_input, cell_x, cell_y, imsizex, imsizey,
213  oversample_ratio);
214  return init_degrid_operator_2d<T>(uv_vis.u, uv_vis.v, uv_vis.w, uv_vis.weights, imsizey, imsizex,
215  oversample_ratio, kernel, Ju, Jw, w_stacking, cell_x, cell_y,
216  absolute_error, relative_error, dde);
217 }

References purify::utilities::convert_to_pixels(), and purify::utilities::w_stacking().

◆ init_degrid_operator_2d() [3/4]

template<class T >
std::shared_ptr<sopt::LinearTransform<T> > purify::measurementoperator::init_degrid_operator_2d ( const Vector< t_real > &  u,
const Vector< t_real > &  v,
const Vector< t_real > &  w,
const Vector< t_complex > &  weights,
const t_uint &  imsizey,
const t_uint &  imsizex,
const t_real &  oversample_ratio = 2,
const kernels::kernel  kernel = kernels::kernel::kb,
const t_uint  Ju = 4,
const t_uint  Jv = 4,
const bool  w_stacking = false,
const t_real &  cellx = 1,
const t_real &  celly = 1 
)

Returns linear transform that is the standard degridding operator.

Definition at line 608 of file operators.h.

613  {
614  const operators::fftw_plan ft_plan = operators::fftw_plan::measure;
615  std::array<t_int, 3> N = {0, 1, static_cast<t_int>(imsizey * imsizex)};
616  std::array<t_int, 3> M = {0, 1, static_cast<t_int>(u.size())};
617  sopt::OperatorFunction<T> directDegrid, indirectDegrid;
618  std::tie(directDegrid, indirectDegrid) = purify::operators::base_degrid_operator_2d<T>(
619  u, v, w, weights, imsizey, imsizex, oversample_ratio, kernel, Ju, Jv, ft_plan, w_stacking,
620  cellx, celly);
621  return std::make_shared<sopt::LinearTransform<T>>(directDegrid, M, indirectDegrid, N);
622 }
const std::vector< t_real > u
data for u coordinate
Definition: operators.cc:18
const std::vector< t_real > v
data for v coordinate
Definition: operators.cc:20
fftw_plan
enum for fftw plans
Definition: operators.h:389

References purify::operators::measure, operators_test::u, operators_test::v, and purify::utilities::w_stacking().

Referenced by degrid_operator_ctor(), dirty_visibilities(), main(), purify::factory::measurement_operator_factory(), padmm(), DegridOperatorFixture::SetUp(), and TEST_CASE().

◆ init_degrid_operator_2d() [4/4]

template<class T >
std::shared_ptr<sopt::LinearTransform<T> > purify::measurementoperator::init_degrid_operator_2d ( const Vector< t_real > &  u,
const Vector< t_real > &  v,
const Vector< t_real > &  w,
const Vector< t_complex > &  weights,
const t_uint  imsizey,
const t_uint  imsizex,
const t_real  oversample_ratio,
const kernels::kernel  kernel,
const t_uint  Ju,
const t_uint  Jw,
const bool  w_stacking,
const t_real  cellx,
const t_real  celly,
const t_real  absolute_error,
const t_real  relative_error,
const dde_type  dde 
)

Returns linear transform that is the standard degridding operator.

Definition at line 189 of file wproj_operators.h.

194  {
195  const operators::fftw_plan ft_plan = operators::fftw_plan::measure;
196  std::array<t_int, 3> N = {0, 1, static_cast<t_int>(imsizey * imsizex)};
197  std::array<t_int, 3> M = {0, 1, static_cast<t_int>(u.size())};
198  sopt::OperatorFunction<T> directDegrid, indirectDegrid;
199  std::tie(directDegrid, indirectDegrid) = purify::operators::base_degrid_operator_2d<T>(
200  u, v, w, weights, imsizey, imsizex, oversample_ratio, kernel, Ju, Jw, ft_plan, w_stacking,
201  cellx, celly, absolute_error, relative_error, dde);
202  auto direct = directDegrid;
203  auto indirect = indirectDegrid;
204  return std::make_shared<sopt::LinearTransform<T>>(direct, M, indirect, N);
205 }

References purify::operators::measure, operators_test::u, operators_test::v, and purify::utilities::w_stacking().