diff --git a/doc/filters.texi b/doc/filters.texi
index 860d1eadca..3c5941d748 100644
--- a/doc/filters.texi
+++ b/doc/filters.texi
@@ -13675,6 +13675,19 @@ Xmap and Ymap input video streams must be of same dimensions. Output video strea
 will have Xmap/Ymap video stream dimensions.
 Xmap and Ymap input video streams are 16bit depth, single channel.
 
+@table @option
+@item interpolation
+
+Select interpolation mode.
+
+@table @samp
+@item nearest
+Use values from the nearest neighbor interpolation.
+@item bilinear
+Interpolate values using the bilinear interpolation.
+@end table
+@end table
+
 @section removegrain
 
 The removegrain filter is a spatial denoiser for progressive video.
diff --git a/libavfilter/vf_remap.c b/libavfilter/vf_remap.c
index 48ec38af7c..9dfb9ce0be 100644
--- a/libavfilter/vf_remap.c
+++ b/libavfilter/vf_remap.c
@@ -47,6 +47,7 @@
 
 typedef struct RemapContext {
     const AVClass *class;
+    int interpolation;
     int nb_planes;
     int nb_components;
     int step;
@@ -59,6 +60,9 @@ typedef struct RemapContext {
 #define FLAGS AV_OPT_FLAG_FILTERING_PARAM|AV_OPT_FLAG_VIDEO_PARAM
 
 static const AVOption remap_options[] = {
+    { "interpolation", "select interpolation mode", OFFSET(interpolation),    AV_OPT_TYPE_INT  , {.i64=0}, 0, 1, FLAGS, "interp" },
+        { "nearest",  "use values from the nearest defined points", 0, AV_OPT_TYPE_CONST, {.i64=0}, 0, 0, FLAGS, "interp" },
+        { "bilinear", "use values from the linear interpolation",   0, AV_OPT_TYPE_CONST, {.i64=1}, 0, 0, FLAGS, "interp" },
     { NULL }
 };
 
@@ -121,6 +125,138 @@ fail:
     return ret;
 }
 
+static av_always_inline float lerpf(float s, float e, float t)
+{
+    return s + (e - s) * t;
+}
+
+static av_always_inline int blerp(int c00, int c10,
+                                  int c01, int c11,
+                                  float tx, float ty)
+{
+    return lerpf(lerpf(c00, c10, tx), lerpf(c01, c11, tx), ty);
+}
+
+#define DEFINE_INTERP_NEAREST_PLANAR(bits)                                              \
+static av_always_inline int interp_nearest_planar##bits(const uint##bits##_t *src,      \
+                                                   int slinesize,                       \
+                                                   const uint16_t *xmap, int xlinesize, \
+                                                   const uint16_t *ymap, int ylinesize, \
+                                                   int w, int h, int x, int y)          \
+{                                                                                       \
+    if (ymap[x] < h && xmap[x] < w) {                                                   \
+        return src[ymap[x] * slinesize + xmap[x]];                                      \
+    }                                                                                   \
+    return 0;                                                                           \
+}
+
+#define DEFINE_INTERP_NEAREST_PACKED(bits)                                              \
+static av_always_inline int interp_nearest_packed##bits(const uint##bits##_t *src,      \
+                                                   int slinesize,                       \
+                                                   const uint16_t *xmap, int xlinesize, \
+                                                   const uint16_t *ymap, int ylinesize, \
+                                                   int w, int h, int x, int y,          \
+                                                   int step, int c)                     \
+{                                                                                       \
+    if (ymap[x] < h && xmap[x] < w) {                                                   \
+        return src[ymap[x] * slinesize + xmap[x] * step + c];                           \
+    }                                                                                   \
+    return 0;                                                                           \
+}
+
+#define DEFINE_INTERP_BILINEAR_PLANAR(bits)                                             \
+static av_always_inline int interp_bilinear_planar##bits(const uint##bits##_t *src,     \
+                                                   int slinesize,                       \
+                                                   const uint16_t *xmap, int xlinesize, \
+                                                   const uint16_t *ymap, int ylinesize, \
+                                                   int w, int h, int x, int y)          \
+{                                                                                       \
+    int c00, c10, c01, c11;                                                             \
+    float xa, ya, ymin, ymax, xmin, xmax;                                               \
+    const int x00 = xmap[x];                                                            \
+    const int x10 = xmap[x + (y < h - 1 ? xlinesize : 0)];                              \
+    const int x01 = xmap[FFMIN(x + 1, w-1)];                                            \
+    const int x11 = xmap[FFMIN(x + 1, w-1) + (y < h - 1 ? xlinesize : 0)];              \
+                                                                                        \
+    const int y00 = ymap[x];                                                            \
+    const int y10 = ymap[x + (y < h - 1 ? ylinesize : 0)];                              \
+    const int y01 = ymap[FFMIN(x + 1, w-1)];                                            \
+    const int y11 = ymap[FFMIN(x + 1, w-1) + (y < h - 1 ? ylinesize : 0)];              \
+                                                                                        \
+    if (x00 >= w || x10 >= w || x01 >= w || x11 >= w ||                                 \
+        y00 >= h || y10 >= h || y01 >= h || y11 >= h)                                   \
+        return 0;                                                                       \
+                                                                                        \
+    c00 = src[y00 * slinesize + x00];                                                   \
+    c10 = src[y10 * slinesize + x10];                                                   \
+    c01 = src[y01 * slinesize + x01];                                                   \
+    c11 = src[y11 * slinesize + x11];                                                   \
+                                                                                        \
+    xa = (x00 + x10 + x01 + x11) / 4.f;                                                 \
+    ya = (y00 + y10 + y01 + y11) / 4.f;                                                 \
+                                                                                        \
+    xmin = FFMIN(FFMIN(x00, x11), FFMIN(x01, x10));                                     \
+    ymin = FFMIN(FFMIN(y00, y11), FFMIN(y01, y10));                                     \
+    xmax = FFMAX(FFMAX(x00, x11), FFMAX(x01, x10));                                     \
+    ymax = FFMAX(FFMAX(y00, y11), FFMAX(y01, y10));                                     \
+                                                                                        \
+    return blerp(c00, c10, c01, c11, (xa - xmin) / FFMAX((xmax - xmin), 1.f),           \
+                                     (ya - ymin) / FFMAX((ymax - ymin), 1.f));          \
+}
+
+#define DEFINE_INTERP_BILINEAR_PACKED(bits)                                             \
+static av_always_inline int interp_bilinear_packed##bits(const uint##bits##_t *src,     \
+                                                   int slinesize,                       \
+                                                   const uint16_t *xmap, int xlinesize, \
+                                                   const uint16_t *ymap, int ylinesize, \
+                                                   int w, int h, int x, int y,          \
+                                                   int step, int c)                     \
+{                                                                                       \
+    int c00, c10, c01, c11;                                                             \
+    float xa, ya, ymin, ymax, xmin, xmax;                                               \
+    const int x00 = xmap[x];                                                            \
+    const int x10 = xmap[x + (y < h - 1 ? xlinesize : 0)];                              \
+    const int x01 = xmap[FFMIN(x + 1, w-1)];                                            \
+    const int x11 = xmap[FFMIN(x + 1, w-1) + (y < h - 1 ? xlinesize : 0)];              \
+                                                                                        \
+    const int y00 = ymap[x];                                                            \
+    const int y10 = ymap[x + (y < h - 1 ? ylinesize : 0)];                              \
+    const int y01 = ymap[FFMIN(x + 1, w-1)];                                            \
+    const int y11 = ymap[FFMIN(x + 1, w-1) + (y < h - 1 ? ylinesize : 0)];              \
+                                                                                        \
+    if (x00 >= w || x10 >= w || x01 >= w || x11 >= w ||                                 \
+        y00 >= h || y10 >= h || y01 >= h || y11 >= h)                                   \
+        return 0;                                                                       \
+                                                                                        \
+    c00 = src[y00 * slinesize + x00 * step + c];                                        \
+    c10 = src[y10 * slinesize + x10 * step + c];                                        \
+    c01 = src[y01 * slinesize + x01 * step + c];                                        \
+    c11 = src[y11 * slinesize + x11 * step + c];                                        \
+                                                                                        \
+    xa = (x00 + x10 + x01 + x11) / 4.f;                                                 \
+    ya = (y00 + y10 + y01 + y11) / 4.f;                                                 \
+                                                                                        \
+    xmin = FFMIN(FFMIN(x00, x11), FFMIN(x01, x10));                                     \
+    ymin = FFMIN(FFMIN(y00, y11), FFMIN(y01, y10));                                     \
+    xmax = FFMAX(FFMAX(x00, x11), FFMAX(x01, x10));                                     \
+    ymax = FFMAX(FFMAX(y00, y11), FFMAX(y01, y10));                                     \
+                                                                                        \
+    return blerp(c00, c10, c01, c11, (xa - xmin) / FFMAX((xmax - xmin), 1.f),           \
+                                     (ya - ymin) / FFMAX((ymax - ymin), 1.f));          \
+}
+
+DEFINE_INTERP_NEAREST_PLANAR(8)
+DEFINE_INTERP_NEAREST_PLANAR(16)
+
+DEFINE_INTERP_BILINEAR_PLANAR(8)
+DEFINE_INTERP_BILINEAR_PLANAR(16)
+
+DEFINE_INTERP_NEAREST_PACKED(8)
+DEFINE_INTERP_NEAREST_PACKED(16)
+
+DEFINE_INTERP_BILINEAR_PACKED(8)
+DEFINE_INTERP_BILINEAR_PACKED(16)
+
 /**
  * remap_planar algorithm expects planes of same size
  * pixels are copied from source to target using :
@@ -151,11 +287,9 @@ static int remap_planar##bits##_##name##_slice(AVFilterContext *ctx, void *arg,
                                                                                             \
         for (y = slice_start; y < slice_end; y++) {                                         \
             for (x = 0; x < out->width; x++) {                                              \
-                if (ymap[x] < in->height && xmap[x] < in->width) {                          \
-                    dst[x] = src[ymap[x] * slinesize + xmap[x]];                            \
-                } else {                                                                    \
-                    dst[x] = 0;                                                             \
-                }                                                                           \
+                dst[x] = interp_##name##_planar##bits(src, slinesize, xmap, xlinesize,      \
+                                                      ymap, ylinesize, in->width,           \
+                                                      in->height, x, y);                    \
             }                                                                               \
             dst  += dlinesize;                                                              \
             xmap += xlinesize;                                                              \
@@ -169,6 +303,9 @@ static int remap_planar##bits##_##name##_slice(AVFilterContext *ctx, void *arg,
 DEFINE_REMAP_PLANAR_FUNC(nearest, 8, 1)
 DEFINE_REMAP_PLANAR_FUNC(nearest, 16, 2)
 
+DEFINE_REMAP_PLANAR_FUNC(bilinear, 8, 1)
+DEFINE_REMAP_PLANAR_FUNC(bilinear, 16, 2)
+
 /**
  * remap_packed algorithm expects pixels with both padded bits (step) and
  * number of components correctly set.
@@ -200,11 +337,10 @@ static int remap_packed##bits##_##name##_slice(AVFilterContext *ctx, void *arg,
     for (y = slice_start; y < slice_end; y++) {                                             \
         for (x = 0; x < out->width; x++) {                                                  \
             for (c = 0; c < td->nb_components; c++) {                                       \
-                if (ymap[x] < in->height && xmap[x] < in->width) {                          \
-                    dst[x * step + c] = src[ymap[x] * slinesize + xmap[x] * step + c];      \
-                } else {                                                                    \
-                    dst[x * step + c] = 0;                                                  \
-                }                                                                           \
+                dst[x * step + c] = interp_##name##_packed##bits(src, slinesize,            \
+                                                      xmap, xlinesize,                      \
+                                                      ymap, ylinesize, in->width,           \
+                                                      in->height, x, y, step, c);           \
             }                                                                               \
         }                                                                                   \
         dst  += dlinesize;                                                                  \
@@ -218,6 +354,9 @@ static int remap_packed##bits##_##name##_slice(AVFilterContext *ctx, void *arg,
 DEFINE_REMAP_PACKED_FUNC(nearest, 8, 1)
 DEFINE_REMAP_PACKED_FUNC(nearest, 16, 2)
 
+DEFINE_REMAP_PACKED_FUNC(bilinear, 8, 1)
+DEFINE_REMAP_PACKED_FUNC(bilinear, 16, 2)
+
 static int config_input(AVFilterLink *inlink)
 {
     AVFilterContext *ctx = inlink->dst;
@@ -229,15 +368,15 @@ static int config_input(AVFilterLink *inlink)
 
     if (desc->comp[0].depth == 8) {
         if (s->nb_planes > 1 || s->nb_components == 1) {
-            s->remap_slice = remap_planar8_nearest_slice;
+            s->remap_slice = s->interpolation ? remap_planar8_bilinear_slice : remap_planar8_nearest_slice;
         } else {
-            s->remap_slice = remap_packed8_nearest_slice;
+            s->remap_slice = s->interpolation ? remap_packed8_bilinear_slice : remap_packed8_nearest_slice;
         }
     } else {
         if (s->nb_planes > 1 || s->nb_components == 1) {
-            s->remap_slice = remap_planar16_nearest_slice;
+            s->remap_slice = s->interpolation ? remap_planar16_bilinear_slice : remap_planar16_nearest_slice;
         } else {
-            s->remap_slice = remap_packed16_nearest_slice;
+            s->remap_slice = s->interpolation ? remap_packed16_bilinear_slice : remap_packed16_nearest_slice;
         }
     }
 
