diff --git a/include/gemmi/grid.hpp b/include/gemmi/grid.hpp
index 1dc2e811..0be6b71e 100644
--- a/include/gemmi/grid.hpp
+++ b/include/gemmi/grid.hpp
@@ -783,20 +783,21 @@ struct Grid : GridBase<T> {
     for (T& x : data)
       x = static_cast<T>((x - stats.dmean) / stats.rms);
   }
-
-  // TODO: can it be replaced with interpolate_grid(dest, src, Transform(), order)?
-  void resample_to(Grid<T>& dest, int order) const {
-    dest.check_not_empty();
-    int idx = 0;
-    for (int w = 0; w < dest.nw; ++w)
-      for (int v = 0; v < dest.nv; ++v)
-        for (int u = 0; u < dest.nu; ++u, ++idx) {
-          const Fractional f = dest.get_fractional(u, v, w);
-          dest.data[idx] = interpolate(f, order);
-        }
-  }
 };
 
+// TODO: add argument Box<Fractional> src_extent
+template<typename T>
+void interpolate_grid(Grid<T>& dest, const Grid<T>& src, const Transform& tr, int order=2) {
+  FTransform frac_tr = src.unit_cell.frac.combine(tr).combine(dest.unit_cell.orth);
+  size_t idx = 0;
+  for (int w = 0; w != dest.nw; ++w)
+    for (int v = 0; v != dest.nv; ++v)
+      for (int u = 0; u != dest.nu; ++u, ++idx) {
+        Fractional dest_fr = dest.get_fractional(u, v, w);
+        Fractional src_fr = frac_tr.apply(dest_fr);
+        dest.data[idx] = src.interpolate(src_fr, order);
+      }
+}
 
 template<typename T>
 Correlation calculate_correlation(const GridBase<T>& a, const GridBase<T>& b) {
diff --git a/include/gemmi/solmask.hpp b/include/gemmi/solmask.hpp
index 6e0eebef..0e76ca71 100644
--- a/include/gemmi/solmask.hpp
+++ b/include/gemmi/solmask.hpp
@@ -332,20 +332,6 @@ struct SolventMasker {
 #endif
 };
 
-// TODO: add argument Box<Fractional> src_extent
-template<typename T>
-void interpolate_grid(Grid<T>& dest, const Grid<T>& src, const Transform& tr, int order=2) {
-  FTransform frac_tr = src.unit_cell.frac.combine(tr).combine(dest.unit_cell.orth);
-  size_t idx = 0;
-  for (int w = 0; w != dest.nw; ++w)
-    for (int v = 0; v != dest.nv; ++v)
-      for (int u = 0; u != dest.nu; ++u, ++idx) {
-        Fractional dest_fr = dest.get_fractional(u, v, w);
-        Fractional src_fr = frac_tr.apply(dest_fr);
-        dest.data[idx] = src.interpolate(src_fr, order);
-      }
-}
-
 struct NodeInfo {
   double dist_sq;  // distance from the nearest atom
   bool found = false;  // the mask flag
diff --git a/python/grid.cpp b/python/grid.cpp
index cd990f9a..e7ff6940 100644
--- a/python/grid.cpp
+++ b/python/grid.cpp
@@ -121,7 +121,6 @@ nb::class_<Grid<T>, GridBase<T>> add_grid_common(nb::module_& m, const std::stri
     .def("symmetrize_max", &Gr::symmetrize_max)
     .def("symmetrize_abs_max", &Gr::symmetrize_abs_max)
     .def("symmetrize_sum", &Gr::symmetrize_sum)
-    .def("resample_to", &Gr::resample_to, nb::arg("dest"), nb::arg("order"))
     .def("masked_asu", &masked_asu<T>, nb::keep_alive<0, 1>())
     .def("mask_points_in_constant_radius", &mask_points_in_constant_radius<T>,
          nb::arg("model"), nb::arg("radius"), nb::arg("value"),