diff --git a/bfps/cpp/interpolator.cpp b/bfps/cpp/interpolator.cpp
index 254213e21ac96545ff46669902d3cde64f3a99a6..4f1467952dfed988f2d854cc3130611c9dadda5c 100644
--- a/bfps/cpp/interpolator.cpp
+++ b/bfps/cpp/interpolator.cpp
@@ -146,8 +146,6 @@ void interpolator<rnumber, interp_neighbours>::operator()(
         this->compute_beta(deriv[2], xx[2], bz);
     }
     std::fill_n(dest, 3, 0);
-    rnumber *ff0 = this->f0 + this->buffer_size;
-    rnumber *ff1 = this->f1 + this->buffer_size;
     for (int iz = -interp_neighbours; iz <= interp_neighbours+1; iz++)
     for (int iy = -interp_neighbours; iy <= interp_neighbours+1; iy++)
     for (int ix = -interp_neighbours; ix <= interp_neighbours+1; ix++)
@@ -155,8 +153,10 @@ void interpolator<rnumber, interp_neighbours>::operator()(
         {
             ptrdiff_t tindex = ((ptrdiff_t(    xg[2]+iz                             ) *this->descriptor->sizes[1] +
                                  ptrdiff_t(MOD(xg[1]+iy, this->descriptor->sizes[1])))*this->descriptor->sizes[2] +
-                                 ptrdiff_t(MOD(xg[0]+ix, this->descriptor->sizes[2])))*3+c;
-            dest[c] += (ff0[tindex]*(1-t) + t*ff1[tindex])*(bz[iz+interp_neighbours]*by[iy+interp_neighbours]*bx[ix+interp_neighbours]);
+                                 ptrdiff_t(MOD(xg[0]+ix, this->descriptor->sizes[2])))*3+c + this->buffer_size;
+            dest[c] += (this->f0[tindex]*(1-t) + t*this->f1[tindex])*(bz[iz+interp_neighbours]*
+                                                                      by[iy+interp_neighbours]*
+                                                                      bx[ix+interp_neighbours]);
         }
 }