git.sesse.net Git - ffmpeg/blob - libavfilter/vf_v360.c

   1 /*
   2  * Copyright (c) 2019 Eugene Lyapustin
   3  *
   4  * This file is part of FFmpeg.
   5  *
   6  * FFmpeg is free software; you can redistribute it and/or
   7  * modify it under the terms of the GNU Lesser General Public
   8  * License as published by the Free Software Foundation; either
   9  * version 2.1 of the License, or (at your option) any later version.
  10  *
  11  * FFmpeg is distributed in the hope that it will be useful,
  12  * but WITHOUT ANY WARRANTY; without even the implied warranty of
  13  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  14  * Lesser General Public License for more details.
  15  *
  16  * You should have received a copy of the GNU Lesser General Public
  17  * License along with FFmpeg; if not, write to the Free Software
  18  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA 02110-1301 USA
  19  */
  20
  21 /**
  22  * @file
  23  * 360 video conversion filter.
  24  * Principle of operation:
  25  *
  26  * (for each pixel in output frame)
  27  * 1) Calculate OpenGL-like coordinates (x, y, z) for pixel position (i, j)
  28  * 2) Apply 360 operations (rotation, mirror) to (x, y, z)
  29  * 3) Calculate pixel position (u, v) in input frame
  30  * 4) Calculate interpolation window and weight for each pixel
  31  *
  32  * (for each frame)
  33  * 5) Remap input frame to output frame using precalculated data
  34  */
  35
  36 #include "libavutil/avassert.h"
  37 #include "libavutil/imgutils.h"
  38 #include "libavutil/pixdesc.h"
  39 #include "libavutil/opt.h"
  40 #include "avfilter.h"
  41 #include "formats.h"
  42 #include "internal.h"
  43 #include "video.h"
  44 #include "v360.h"
  45
  46 typedef struct ThreadData {
  47     AVFrame *in;
  48     AVFrame *out;
  49 } ThreadData;
  50
  51 #define OFFSET(x) offsetof(V360Context, x)
  52 #define FLAGS AV_OPT_FLAG_FILTERING_PARAM|AV_OPT_FLAG_VIDEO_PARAM
  53
  54 static const AVOption v360_options[] = {
  55     {     "input", "set input projection",              OFFSET(in), AV_OPT_TYPE_INT,    {.i64=EQUIRECTANGULAR}, 0,    NB_PROJECTIONS-1, FLAGS, "in" },
  56     {         "e", "equirectangular",                            0, AV_OPT_TYPE_CONST,  {.i64=EQUIRECTANGULAR}, 0,                   0, FLAGS, "in" },
  57     {  "equirect", "equirectangular",                            0, AV_OPT_TYPE_CONST,  {.i64=EQUIRECTANGULAR}, 0,                   0, FLAGS, "in" },
  58     {      "c3x2", "cubemap 3x2",                                0, AV_OPT_TYPE_CONST,  {.i64=CUBEMAP_3_2},     0,                   0, FLAGS, "in" },
  59     {      "c6x1", "cubemap 6x1",                                0, AV_OPT_TYPE_CONST,  {.i64=CUBEMAP_6_1},     0,                   0, FLAGS, "in" },
  60     {       "eac", "equi-angular cubemap",                       0, AV_OPT_TYPE_CONST,  {.i64=EQUIANGULAR},     0,                   0, FLAGS, "in" },
  61     {  "dfisheye", "dual fisheye",                               0, AV_OPT_TYPE_CONST,  {.i64=DUAL_FISHEYE},    0,                   0, FLAGS, "in" },
  62     {    "barrel", "barrel facebook's 360 format",               0, AV_OPT_TYPE_CONST,  {.i64=BARREL},          0,                   0, FLAGS, "in" },
  63     {        "fb", "barrel facebook's 360 format",               0, AV_OPT_TYPE_CONST,  {.i64=BARREL},          0,                   0, FLAGS, "in" },
  64     {      "c1x6", "cubemap 1x6",                                0, AV_OPT_TYPE_CONST,  {.i64=CUBEMAP_1_6},     0,                   0, FLAGS, "in" },
  65     {    "output", "set output projection",            OFFSET(out), AV_OPT_TYPE_INT,    {.i64=CUBEMAP_3_2},     0,    NB_PROJECTIONS-1, FLAGS, "out" },
  66     {         "e", "equirectangular",                            0, AV_OPT_TYPE_CONST,  {.i64=EQUIRECTANGULAR}, 0,                   0, FLAGS, "out" },
  67     {  "equirect", "equirectangular",                            0, AV_OPT_TYPE_CONST,  {.i64=EQUIRECTANGULAR}, 0,                   0, FLAGS, "out" },
  68     {      "c3x2", "cubemap 3x2",                                0, AV_OPT_TYPE_CONST,  {.i64=CUBEMAP_3_2},     0,                   0, FLAGS, "out" },
  69     {      "c6x1", "cubemap 6x1",                                0, AV_OPT_TYPE_CONST,  {.i64=CUBEMAP_6_1},     0,                   0, FLAGS, "out" },
  70     {       "eac", "equi-angular cubemap",                       0, AV_OPT_TYPE_CONST,  {.i64=EQUIANGULAR},     0,                   0, FLAGS, "out" },
  71     {      "flat", "regular video",                              0, AV_OPT_TYPE_CONST,  {.i64=FLAT},            0,                   0, FLAGS, "out" },
  72     {"rectilinear", "regular video",                             0, AV_OPT_TYPE_CONST,  {.i64=FLAT},            0,                   0, FLAGS, "out" },
  73     {  "gnomonic", "regular video",                              0, AV_OPT_TYPE_CONST,  {.i64=FLAT},            0,                   0, FLAGS, "out" },
  74     {    "barrel", "barrel facebook's 360 format",               0, AV_OPT_TYPE_CONST,  {.i64=BARREL},          0,                   0, FLAGS, "out" },
  75     {        "fb", "barrel facebook's 360 format",               0, AV_OPT_TYPE_CONST,  {.i64=BARREL},          0,                   0, FLAGS, "out" },
  76     {      "c1x6", "cubemap 1x6",                                0, AV_OPT_TYPE_CONST,  {.i64=CUBEMAP_1_6},     0,                   0, FLAGS, "out" },
  77     {        "sg", "stereographic",                              0, AV_OPT_TYPE_CONST,  {.i64=STEREOGRAPHIC},   0,                   0, FLAGS, "out" },
  78     {    "interp", "set interpolation method",      OFFSET(interp), AV_OPT_TYPE_INT,    {.i64=BILINEAR},        0, NB_INTERP_METHODS-1, FLAGS, "interp" },
  79     {      "near", "nearest neighbour",                          0, AV_OPT_TYPE_CONST,  {.i64=NEAREST},         0,                   0, FLAGS, "interp" },
  80     {   "nearest", "nearest neighbour",                          0, AV_OPT_TYPE_CONST,  {.i64=NEAREST},         0,                   0, FLAGS, "interp" },
  81     {      "line", "bilinear interpolation",                     0, AV_OPT_TYPE_CONST,  {.i64=BILINEAR},        0,                   0, FLAGS, "interp" },
  82     {    "linear", "bilinear interpolation",                     0, AV_OPT_TYPE_CONST,  {.i64=BILINEAR},        0,                   0, FLAGS, "interp" },
  83     {      "cube", "bicubic interpolation",                      0, AV_OPT_TYPE_CONST,  {.i64=BICUBIC},         0,                   0, FLAGS, "interp" },
  84     {     "cubic", "bicubic interpolation",                      0, AV_OPT_TYPE_CONST,  {.i64=BICUBIC},         0,                   0, FLAGS, "interp" },
  85     {      "lanc", "lanczos interpolation",                      0, AV_OPT_TYPE_CONST,  {.i64=LANCZOS},         0,                   0, FLAGS, "interp" },
  86     {   "lanczos", "lanczos interpolation",                      0, AV_OPT_TYPE_CONST,  {.i64=LANCZOS},         0,                   0, FLAGS, "interp" },
  87     {         "w", "output width",                   OFFSET(width), AV_OPT_TYPE_INT,    {.i64=0},               0,           INT16_MAX, FLAGS, "w"},
  88     {         "h", "output height",                 OFFSET(height), AV_OPT_TYPE_INT,    {.i64=0},               0,           INT16_MAX, FLAGS, "h"},
  89     { "in_forder", "input cubemap face order",   OFFSET(in_forder), AV_OPT_TYPE_STRING, {.str="rludfb"},        0,     NB_DIRECTIONS-1, FLAGS, "in_forder"},
  90     {"out_forder", "output cubemap face order", OFFSET(out_forder), AV_OPT_TYPE_STRING, {.str="rludfb"},        0,     NB_DIRECTIONS-1, FLAGS, "out_forder"},
  91     {   "in_frot", "input cubemap face rotation",  OFFSET(in_frot), AV_OPT_TYPE_STRING, {.str="000000"},        0,     NB_DIRECTIONS-1, FLAGS, "in_frot"},
  92     {  "out_frot", "output cubemap face rotation",OFFSET(out_frot), AV_OPT_TYPE_STRING, {.str="000000"},        0,     NB_DIRECTIONS-1, FLAGS, "out_frot"},
  93     {    "in_pad", "input cubemap pads",            OFFSET(in_pad), AV_OPT_TYPE_FLOAT,  {.dbl=0.f},           0.f,                 1.f, FLAGS, "in_pad"},
  94     {   "out_pad", "output cubemap pads",          OFFSET(out_pad), AV_OPT_TYPE_FLOAT,  {.dbl=0.f},           0.f,                 1.f, FLAGS, "out_pad"},
  95     {       "yaw", "yaw rotation",                     OFFSET(yaw), AV_OPT_TYPE_FLOAT,  {.dbl=0.f},        -180.f,               180.f, FLAGS, "yaw"},
  96     {     "pitch", "pitch rotation",                 OFFSET(pitch), AV_OPT_TYPE_FLOAT,  {.dbl=0.f},        -180.f,               180.f, FLAGS, "pitch"},
  97     {      "roll", "roll rotation",                   OFFSET(roll), AV_OPT_TYPE_FLOAT,  {.dbl=0.f},        -180.f,               180.f, FLAGS, "roll"},
  98     {    "rorder", "rotation order",                OFFSET(rorder), AV_OPT_TYPE_STRING, {.str="ypr"},           0,                   0, FLAGS, "rorder"},
  99     {     "h_fov", "horizontal field of view",       OFFSET(h_fov), AV_OPT_TYPE_FLOAT,  {.dbl=90.f},     0.00001f,               180.f, FLAGS, "h_fov"},
 100     {     "v_fov", "vertical field of view",         OFFSET(v_fov), AV_OPT_TYPE_FLOAT,  {.dbl=45.f},     0.00001f,                90.f, FLAGS, "v_fov"},
 101     {    "h_flip", "flip out video horizontally",   OFFSET(h_flip), AV_OPT_TYPE_BOOL,   {.i64=0},               0,                   1, FLAGS, "h_flip"},
 102     {    "v_flip", "flip out video vertically",     OFFSET(v_flip), AV_OPT_TYPE_BOOL,   {.i64=0},               0,                   1, FLAGS, "v_flip"},
 103     {    "d_flip", "flip out video indepth",        OFFSET(d_flip), AV_OPT_TYPE_BOOL,   {.i64=0},               0,                   1, FLAGS, "d_flip"},
 104     {   "ih_flip", "flip in video horizontally",   OFFSET(ih_flip), AV_OPT_TYPE_BOOL,   {.i64=0},               0,                   1, FLAGS, "ih_flip"},
 105     {   "iv_flip", "flip in video vertically",     OFFSET(iv_flip), AV_OPT_TYPE_BOOL,   {.i64=0},               0,                   1, FLAGS, "iv_flip"},
 106     {  "in_trans", "transpose video input",   OFFSET(in_transpose), AV_OPT_TYPE_BOOL,   {.i64=0},               0,                   1, FLAGS, "in_transpose"},
 107     { "out_trans", "transpose video output", OFFSET(out_transpose), AV_OPT_TYPE_BOOL,   {.i64=0},               0,                   1, FLAGS, "out_transpose"},
 108     { NULL }
 109 };
 110
 111 AVFILTER_DEFINE_CLASS(v360);
 112
 113 static int query_formats(AVFilterContext *ctx)
 114 {
 115     static const enum AVPixelFormat pix_fmts[] = {
 116         // YUVA444
 117         AV_PIX_FMT_YUVA444P,   AV_PIX_FMT_YUVA444P9,
 118         AV_PIX_FMT_YUVA444P10, AV_PIX_FMT_YUVA444P12,
 119         AV_PIX_FMT_YUVA444P16,
 120
 121         // YUVA422
 122         AV_PIX_FMT_YUVA422P,   AV_PIX_FMT_YUVA422P9,
 123         AV_PIX_FMT_YUVA422P10, AV_PIX_FMT_YUVA422P12,
 124         AV_PIX_FMT_YUVA422P16,
 125
 126         // YUVA420
 127         AV_PIX_FMT_YUVA420P,   AV_PIX_FMT_YUVA420P9,
 128         AV_PIX_FMT_YUVA420P10, AV_PIX_FMT_YUVA420P16,
 129
 130         // YUVJ
 131         AV_PIX_FMT_YUVJ444P, AV_PIX_FMT_YUVJ440P,
 132         AV_PIX_FMT_YUVJ422P, AV_PIX_FMT_YUVJ420P,
 133         AV_PIX_FMT_YUVJ411P,
 134
 135         // YUV444
 136         AV_PIX_FMT_YUV444P,   AV_PIX_FMT_YUV444P9,
 137         AV_PIX_FMT_YUV444P10, AV_PIX_FMT_YUV444P12,
 138         AV_PIX_FMT_YUV444P14, AV_PIX_FMT_YUV444P16,
 139
 140         // YUV440
 141         AV_PIX_FMT_YUV440P, AV_PIX_FMT_YUV440P10,
 142         AV_PIX_FMT_YUV440P12,
 143
 144         // YUV422
 145         AV_PIX_FMT_YUV422P,   AV_PIX_FMT_YUV422P9,
 146         AV_PIX_FMT_YUV422P10, AV_PIX_FMT_YUV422P12,
 147         AV_PIX_FMT_YUV422P14, AV_PIX_FMT_YUV422P16,
 148
 149         // YUV420
 150         AV_PIX_FMT_YUV420P,   AV_PIX_FMT_YUV420P9,
 151         AV_PIX_FMT_YUV420P10, AV_PIX_FMT_YUV420P12,
 152         AV_PIX_FMT_YUV420P14, AV_PIX_FMT_YUV420P16,
 153
 154         // YUV411
 155         AV_PIX_FMT_YUV411P,
 156
 157         // YUV410
 158         AV_PIX_FMT_YUV410P,
 159
 160         // GBR
 161         AV_PIX_FMT_GBRP,   AV_PIX_FMT_GBRP9,
 162         AV_PIX_FMT_GBRP10, AV_PIX_FMT_GBRP12,
 163         AV_PIX_FMT_GBRP14, AV_PIX_FMT_GBRP16,
 164
 165         // GBRA
 166         AV_PIX_FMT_GBRAP,   AV_PIX_FMT_GBRAP10,
 167         AV_PIX_FMT_GBRAP12, AV_PIX_FMT_GBRAP16,
 168
 169         // GRAY
 170         AV_PIX_FMT_GRAY8,  AV_PIX_FMT_GRAY9,
 171         AV_PIX_FMT_GRAY10, AV_PIX_FMT_GRAY12,
 172         AV_PIX_FMT_GRAY14, AV_PIX_FMT_GRAY16,
 173
 174         AV_PIX_FMT_NONE
 175     };
 176
 177     AVFilterFormats *fmts_list = ff_make_format_list(pix_fmts);
 178     if (!fmts_list)
 179         return AVERROR(ENOMEM);
 180     return ff_set_common_formats(ctx, fmts_list);
 181 }
 182
 183 #define DEFINE_REMAP1_LINE(bits, div)                                                           \
 184 static void remap1_##bits##bit_line_c(uint8_t *dst, int width, const uint8_t *src,              \
 185                                       ptrdiff_t in_linesize,                                    \
 186                                       const uint16_t *u, const uint16_t *v, const int16_t *ker) \
 187 {                                                                                               \
 188     const uint##bits##_t *s = (const uint##bits##_t *)src;                                      \
 189     uint##bits##_t *d = (uint##bits##_t *)dst;                                                  \
 190                                                                                                 \
 191     in_linesize /= div;                                                                         \
 192                                                                                                 \
 193     for (int x = 0; x < width; x++)                                                             \
 194         d[x] = s[v[x] * in_linesize + u[x]];                                                    \
 195 }
 196
 197 DEFINE_REMAP1_LINE( 8, 1)
 198 DEFINE_REMAP1_LINE(16, 2)
 199
 200 typedef struct XYRemap {
 201     uint16_t u[4][4];
 202     uint16_t v[4][4];
 203     float ker[4][4];
 204 } XYRemap;
 205
 206 /**
 207  * Generate remapping function with a given window size and pixel depth.
 208  *
 209  * @param ws size of interpolation window
 210  * @param bits number of bits per pixel
 211  */
 212 #define DEFINE_REMAP(ws, bits)                                                                             \
 213 static int remap##ws##_##bits##bit_slice(AVFilterContext *ctx, void *arg, int jobnr, int nb_jobs)          \
 214 {                                                                                                          \
 215     ThreadData *td = (ThreadData*)arg;                                                                     \
 216     const V360Context *s = ctx->priv;                                                                      \
 217     const AVFrame *in = td->in;                                                                            \
 218     AVFrame *out = td->out;                                                                                \
 219                                                                                                            \
 220     for (int plane = 0; plane < s->nb_planes; plane++) {                                                   \
 221         const int in_linesize  = in->linesize[plane];                                                      \
 222         const int out_linesize = out->linesize[plane];                                                     \
 223         const uint8_t *src = in->data[plane];                                                              \
 224         uint8_t *dst = out->data[plane];                                                                   \
 225         const int width = s->planewidth[plane];                                                            \
 226         const int height = s->planeheight[plane];                                                          \
 227                                                                                                            \
 228         const int slice_start = (height *  jobnr     ) / nb_jobs;                                          \
 229         const int slice_end   = (height * (jobnr + 1)) / nb_jobs;                                          \
 230                                                                                                            \
 231         for (int y = slice_start; y < slice_end; y++) {                                                    \
 232             const unsigned map = s->map[plane];                                                            \
 233             const uint16_t *u = s->u[map] + y * width * ws * ws;                                           \
 234             const uint16_t *v = s->v[map] + y * width * ws * ws;                                           \
 235             const int16_t *ker = s->ker[map] + y * width * ws * ws;                                        \
 236                                                                                                            \
 237             s->remap_line(dst + y * out_linesize, width, src, in_linesize, u, v, ker);                     \
 238         }                                                                                                  \
 239     }                                                                                                      \
 240                                                                                                            \
 241     return 0;                                                                                              \
 242 }
 243
 244 DEFINE_REMAP(1,  8)
 245 DEFINE_REMAP(2,  8)
 246 DEFINE_REMAP(4,  8)
 247 DEFINE_REMAP(1, 16)
 248 DEFINE_REMAP(2, 16)
 249 DEFINE_REMAP(4, 16)
 250
 251 #define DEFINE_REMAP_LINE(ws, bits, div)                                                                   \
 252 static void remap##ws##_##bits##bit_line_c(uint8_t *dst, int width, const uint8_t *src,                    \
 253                                            ptrdiff_t in_linesize,                                          \
 254                                            const uint16_t *u, const uint16_t *v, const int16_t *ker)       \
 255 {                                                                                                          \
 256     const uint##bits##_t *s = (const uint##bits##_t *)src;                                                 \
 257     uint##bits##_t *d = (uint##bits##_t *)dst;                                                             \
 258                                                                                                            \
 259     in_linesize /= div;                                                                                    \
 260                                                                                                            \
 261     for (int x = 0; x < width; x++) {                                                                      \
 262         const uint16_t *uu = u + x * ws * ws;                                                              \
 263         const uint16_t *vv = v + x * ws * ws;                                                              \
 264         const int16_t *kker = ker + x * ws * ws;                                                           \
 265         int tmp = 0;                                                                                       \
 266                                                                                                            \
 267         for (int i = 0; i < ws; i++) {                                                                     \
 268             for (int j = 0; j < ws; j++) {                                                                 \
 269                 tmp += kker[i * ws + j] * s[vv[i * ws + j] * in_linesize + uu[i * ws + j]];                \
 270             }                                                                                              \
 271         }                                                                                                  \
 272                                                                                                            \
 273         d[x] = av_clip_uint##bits(tmp >> 14);                                                              \
 274     }                                                                                                      \
 275 }
 276
 277 DEFINE_REMAP_LINE(2,  8, 1)
 278 DEFINE_REMAP_LINE(4,  8, 1)
 279 DEFINE_REMAP_LINE(2, 16, 2)
 280 DEFINE_REMAP_LINE(4, 16, 2)
 281
 282 void ff_v360_init(V360Context *s, int depth)
 283 {
 284     switch (s->interp) {
 285     case NEAREST:
 286         s->remap_line = depth <= 8 ? remap1_8bit_line_c : remap1_16bit_line_c;
 287         break;
 288     case BILINEAR:
 289         s->remap_line = depth <= 8 ? remap2_8bit_line_c : remap2_16bit_line_c;
 290         break;
 291     case BICUBIC:
 292     case LANCZOS:
 293         s->remap_line = depth <= 8 ? remap4_8bit_line_c : remap4_16bit_line_c;
 294         break;
 295     }
 296
 297     if (ARCH_X86)
 298         ff_v360_init_x86(s, depth);
 299 }
 300
 301 /**
 302  * Save nearest pixel coordinates for remapping.
 303  *
 304  * @param du horizontal relative coordinate
 305  * @param dv vertical relative coordinate
 306  * @param r_tmp calculated 4x4 window
 307  * @param u u remap data
 308  * @param v v remap data
 309  * @param ker ker remap data
 310  */
 311 static void nearest_kernel(float du, float dv, const XYRemap *r_tmp,
 312                            uint16_t *u, uint16_t *v, int16_t *ker)
 313 {
 314     const int i = roundf(dv) + 1;
 315     const int j = roundf(du) + 1;
 316
 317     u[0] = r_tmp->u[i][j];
 318     v[0] = r_tmp->v[i][j];
 319 }
 320
 321 /**
 322  * Calculate kernel for bilinear interpolation.
 323  *
 324  * @param du horizontal relative coordinate
 325  * @param dv vertical relative coordinate
 326  * @param r_tmp calculated 4x4 window
 327  * @param u u remap data
 328  * @param v v remap data
 329  * @param ker ker remap data
 330  */
 331 static void bilinear_kernel(float du, float dv, const XYRemap *r_tmp,
 332                             uint16_t *u, uint16_t *v, int16_t *ker)
 333 {
 334     int i, j;
 335
 336     for (i = 0; i < 2; i++) {
 337         for (j = 0; j < 2; j++) {
 338             u[i * 2 + j] = r_tmp->u[i + 1][j + 1];
 339             v[i * 2 + j] = r_tmp->v[i + 1][j + 1];
 340         }
 341     }
 342
 343     ker[0] = (1.f - du) * (1.f - dv) * 16384;
 344     ker[1] =        du  * (1.f - dv) * 16384;
 345     ker[2] = (1.f - du) *        dv  * 16384;
 346     ker[3] =        du  *        dv  * 16384;
 347 }
 348
 349 /**
 350  * Calculate 1-dimensional cubic coefficients.
 351  *
 352  * @param t relative coordinate
 353  * @param coeffs coefficients
 354  */
 355 static inline void calculate_bicubic_coeffs(float t, float *coeffs)
 356 {
 357     const float tt  = t * t;
 358     const float ttt = t * t * t;
 359
 360     coeffs[0] =     - t / 3.f + tt / 2.f - ttt / 6.f;
 361     coeffs[1] = 1.f - t / 2.f - tt       + ttt / 2.f;
 362     coeffs[2] =       t       + tt / 2.f - ttt / 2.f;
 363     coeffs[3] =     - t / 6.f            + ttt / 6.f;
 364 }
 365
 366 /**
 367  * Calculate kernel for bicubic interpolation.
 368  *
 369  * @param du horizontal relative coordinate
 370  * @param dv vertical relative coordinate
 371  * @param r_tmp calculated 4x4 window
 372  * @param u u remap data
 373  * @param v v remap data
 374  * @param ker ker remap data
 375  */
 376 static void bicubic_kernel(float du, float dv, const XYRemap *r_tmp,
 377                            uint16_t *u, uint16_t *v, int16_t *ker)
 378 {
 379     int i, j;
 380     float du_coeffs[4];
 381     float dv_coeffs[4];
 382
 383     calculate_bicubic_coeffs(du, du_coeffs);
 384     calculate_bicubic_coeffs(dv, dv_coeffs);
 385
 386     for (i = 0; i < 4; i++) {
 387         for (j = 0; j < 4; j++) {
 388             u[i * 4 + j] = r_tmp->u[i][j];
 389             v[i * 4 + j] = r_tmp->v[i][j];
 390             ker[i * 4 + j] = du_coeffs[j] * dv_coeffs[i] * 16384;
 391         }
 392     }
 393 }
 394
 395 /**
 396  * Calculate 1-dimensional lanczos coefficients.
 397  *
 398  * @param t relative coordinate
 399  * @param coeffs coefficients
 400  */
 401 static inline void calculate_lanczos_coeffs(float t, float *coeffs)
 402 {
 403     int i;
 404     float sum = 0.f;
 405
 406     for (i = 0; i < 4; i++) {
 407         const float x = M_PI * (t - i + 1);
 408         if (x == 0.f) {
 409             coeffs[i] = 1.f;
 410         } else {
 411             coeffs[i] = sinf(x) * sinf(x / 2.f) / (x * x / 2.f);
 412         }
 413         sum += coeffs[i];
 414     }
 415
 416     for (i = 0; i < 4; i++) {
 417         coeffs[i] /= sum;
 418     }
 419 }
 420
 421 /**
 422  * Calculate kernel for lanczos interpolation.
 423  *
 424  * @param du horizontal relative coordinate
 425  * @param dv vertical relative coordinate
 426  * @param r_tmp calculated 4x4 window
 427  * @param u u remap data
 428  * @param v v remap data
 429  * @param ker ker remap data
 430  */
 431 static void lanczos_kernel(float du, float dv, const XYRemap *r_tmp,
 432                            uint16_t *u, uint16_t *v, int16_t *ker)
 433 {
 434     int i, j;
 435     float du_coeffs[4];
 436     float dv_coeffs[4];
 437
 438     calculate_lanczos_coeffs(du, du_coeffs);
 439     calculate_lanczos_coeffs(dv, dv_coeffs);
 440
 441     for (i = 0; i < 4; i++) {
 442         for (j = 0; j < 4; j++) {
 443             u[i * 4 + j] = r_tmp->u[i][j];
 444             v[i * 4 + j] = r_tmp->v[i][j];
 445             ker[i * 4 + j] = du_coeffs[j] * dv_coeffs[i] * 16384;
 446         }
 447     }
 448 }
 449
 450 /**
 451  * Modulo operation with only positive remainders.
 452  *
 453  * @param a dividend
 454  * @param b divisor
 455  *
 456  * @return positive remainder of (a / b)
 457  */
 458 static inline int mod(int a, int b)
 459 {
 460     const int res = a % b;
 461     if (res < 0) {
 462         return res + b;
 463     } else {
 464         return res;
 465     }
 466 }
 467
 468 /**
 469  * Convert char to corresponding direction.
 470  * Used for cubemap options.
 471  */
 472 static int get_direction(char c)
 473 {
 474     switch (c) {
 475     case 'r':
 476         return RIGHT;
 477     case 'l':
 478         return LEFT;
 479     case 'u':
 480         return UP;
 481     case 'd':
 482         return DOWN;
 483     case 'f':
 484         return FRONT;
 485     case 'b':
 486         return BACK;
 487     default:
 488         return -1;
 489     }
 490 }
 491
 492 /**
 493  * Convert char to corresponding rotation angle.
 494  * Used for cubemap options.
 495  */
 496 static int get_rotation(char c)
 497 {
 498     switch (c) {
 499     case '0':
 500         return ROT_0;
 501     case '1':
 502         return ROT_90;
 503     case '2':
 504         return ROT_180;
 505     case '3':
 506         return ROT_270;
 507     default:
 508         return -1;
 509     }
 510 }
 511
 512 /**
 513  * Convert char to corresponding rotation order.
 514  */
 515 static int get_rorder(char c)
 516 {
 517     switch (c) {
 518     case 'Y':
 519     case 'y':
 520         return YAW;
 521     case 'P':
 522     case 'p':
 523         return PITCH;
 524     case 'R':
 525     case 'r':
 526         return ROLL;
 527     default:
 528         return -1;
 529     }
 530 }
 531
 532 /**
 533  * Prepare data for processing cubemap input format.
 534  *
 535  * @param ctx filter context
 536  *
 537  * @return error code
 538  */
 539 static int prepare_cube_in(AVFilterContext *ctx)
 540 {
 541     V360Context *s = ctx->priv;
 542
 543     for (int face = 0; face < NB_FACES; face++) {
 544         const char c = s->in_forder[face];
 545         int direction;
 546
 547         if (c == '\0') {
 548             av_log(ctx, AV_LOG_ERROR,
 549                    "Incomplete in_forder option. Direction for all 6 faces should be specified.\n");
 550             return AVERROR(EINVAL);
 551         }
 552
 553         direction = get_direction(c);
 554         if (direction == -1) {
 555             av_log(ctx, AV_LOG_ERROR,
 556                    "Incorrect direction symbol '%c' in in_forder option.\n", c);
 557             return AVERROR(EINVAL);
 558         }
 559
 560         s->in_cubemap_face_order[direction] = face;
 561     }
 562
 563     for (int face = 0; face < NB_FACES; face++) {
 564         const char c = s->in_frot[face];
 565         int rotation;
 566
 567         if (c == '\0') {
 568             av_log(ctx, AV_LOG_ERROR,
 569                    "Incomplete in_frot option. Rotation for all 6 faces should be specified.\n");
 570             return AVERROR(EINVAL);
 571         }
 572
 573         rotation = get_rotation(c);
 574         if (rotation == -1) {
 575             av_log(ctx, AV_LOG_ERROR,
 576                    "Incorrect rotation symbol '%c' in in_frot option.\n", c);
 577             return AVERROR(EINVAL);
 578         }
 579
 580         s->in_cubemap_face_rotation[face] = rotation;
 581     }
 582
 583     return 0;
 584 }
 585
 586 /**
 587  * Prepare data for processing cubemap output format.
 588  *
 589  * @param ctx filter context
 590  *
 591  * @return error code
 592  */
 593 static int prepare_cube_out(AVFilterContext *ctx)
 594 {
 595     V360Context *s = ctx->priv;
 596
 597     for (int face = 0; face < NB_FACES; face++) {
 598         const char c = s->out_forder[face];
 599         int direction;
 600
 601         if (c == '\0') {
 602             av_log(ctx, AV_LOG_ERROR,
 603                    "Incomplete out_forder option. Direction for all 6 faces should be specified.\n");
 604             return AVERROR(EINVAL);
 605         }
 606
 607         direction = get_direction(c);
 608         if (direction == -1) {
 609             av_log(ctx, AV_LOG_ERROR,
 610                    "Incorrect direction symbol '%c' in out_forder option.\n", c);
 611             return AVERROR(EINVAL);
 612         }
 613
 614         s->out_cubemap_direction_order[face] = direction;
 615     }
 616
 617     for (int face = 0; face < NB_FACES; face++) {
 618         const char c = s->out_frot[face];
 619         int rotation;
 620
 621         if (c == '\0') {
 622             av_log(ctx, AV_LOG_ERROR,
 623                    "Incomplete out_frot option. Rotation for all 6 faces should be specified.\n");
 624             return AVERROR(EINVAL);
 625         }
 626
 627         rotation = get_rotation(c);
 628         if (rotation == -1) {
 629             av_log(ctx, AV_LOG_ERROR,
 630                    "Incorrect rotation symbol '%c' in out_frot option.\n", c);
 631             return AVERROR(EINVAL);
 632         }
 633
 634         s->out_cubemap_face_rotation[face] = rotation;
 635     }
 636
 637     return 0;
 638 }
 639
 640 static inline void rotate_cube_face(float *uf, float *vf, int rotation)
 641 {
 642     float tmp;
 643
 644     switch (rotation) {
 645     case ROT_0:
 646         break;
 647     case ROT_90:
 648         tmp =  *uf;
 649         *uf = -*vf;
 650         *vf =  tmp;
 651         break;
 652     case ROT_180:
 653         *uf = -*uf;
 654         *vf = -*vf;
 655         break;
 656     case ROT_270:
 657         tmp = -*uf;
 658         *uf =  *vf;
 659         *vf =  tmp;
 660         break;
 661     default:
 662         av_assert0(0);
 663     }
 664 }
 665
 666 static inline void rotate_cube_face_inverse(float *uf, float *vf, int rotation)
 667 {
 668     float tmp;
 669
 670     switch (rotation) {
 671     case ROT_0:
 672         break;
 673     case ROT_90:
 674         tmp = -*uf;
 675         *uf =  *vf;
 676         *vf =  tmp;
 677         break;
 678     case ROT_180:
 679         *uf = -*uf;
 680         *vf = -*vf;
 681         break;
 682     case ROT_270:
 683         tmp =  *uf;
 684         *uf = -*vf;
 685         *vf =  tmp;
 686         break;
 687     default:
 688         av_assert0(0);
 689     }
 690 }
 691
 692 /**
 693  * Normalize vector.
 694  *
 695  * @param vec vector
 696  */
 697 static void normalize_vector(float *vec)
 698 {
 699     const float norm = sqrtf(vec[0] * vec[0] + vec[1] * vec[1] + vec[2] * vec[2]);
 700
 701     vec[0] /= norm;
 702     vec[1] /= norm;
 703     vec[2] /= norm;
 704 }
 705
 706 /**
 707  * Calculate 3D coordinates on sphere for corresponding cubemap position.
 708  * Common operation for every cubemap.
 709  *
 710  * @param s filter context
 711  * @param uf horizontal cubemap coordinate [0, 1)
 712  * @param vf vertical cubemap coordinate [0, 1)
 713  * @param face face of cubemap
 714  * @param vec coordinates on sphere
 715  */
 716 static void cube_to_xyz(const V360Context *s,
 717                         float uf, float vf, int face,
 718                         float *vec)
 719 {
 720     const int direction = s->out_cubemap_direction_order[face];
 721     float l_x, l_y, l_z;
 722
 723     uf /= (1.f - s->out_pad);
 724     vf /= (1.f - s->out_pad);
 725
 726     rotate_cube_face_inverse(&uf, &vf, s->out_cubemap_face_rotation[face]);
 727
 728     switch (direction) {
 729     case RIGHT:
 730         l_x =  1.f;
 731         l_y = -vf;
 732         l_z =  uf;
 733         break;
 734     case LEFT:
 735         l_x = -1.f;
 736         l_y = -vf;
 737         l_z = -uf;
 738         break;
 739     case UP:
 740         l_x =  uf;
 741         l_y =  1.f;
 742         l_z = -vf;
 743         break;
 744     case DOWN:
 745         l_x =  uf;
 746         l_y = -1.f;
 747         l_z =  vf;
 748         break;
 749     case FRONT:
 750         l_x =  uf;
 751         l_y = -vf;
 752         l_z = -1.f;
 753         break;
 754     case BACK:
 755         l_x = -uf;
 756         l_y = -vf;
 757         l_z =  1.f;
 758         break;
 759     }
 760
 761     vec[0] = l_x;
 762     vec[1] = l_y;
 763     vec[2] = l_z;
 764
 765     normalize_vector(vec);
 766 }
 767
 768 /**
 769  * Calculate cubemap position for corresponding 3D coordinates on sphere.
 770  * Common operation for every cubemap.
 771  *
 772  * @param s filter context
 773  * @param vec coordinated on sphere
 774  * @param uf horizontal cubemap coordinate [0, 1)
 775  * @param vf vertical cubemap coordinate [0, 1)
 776  * @param direction direction of view
 777  */
 778 static void xyz_to_cube(const V360Context *s,
 779                         const float *vec,
 780                         float *uf, float *vf, int *direction)
 781 {
 782     const float phi   = atan2f(vec[0], -vec[2]);
 783     const float theta = asinf(-vec[1]);
 784     float phi_norm, theta_threshold;
 785     int face;
 786
 787     if (phi >= -M_PI_4 && phi < M_PI_4) {
 788         *direction = FRONT;
 789         phi_norm = phi;
 790     } else if (phi >= -(M_PI_2 + M_PI_4) && phi < -M_PI_4) {
 791         *direction = LEFT;
 792         phi_norm = phi + M_PI_2;
 793     } else if (phi >= M_PI_4 && phi < M_PI_2 + M_PI_4) {
 794         *direction = RIGHT;
 795         phi_norm = phi - M_PI_2;
 796     } else {
 797         *direction = BACK;
 798         phi_norm = phi + ((phi > 0.f) ? -M_PI : M_PI);
 799     }
 800
 801     theta_threshold = atanf(cosf(phi_norm));
 802     if (theta > theta_threshold) {
 803         *direction = DOWN;
 804     } else if (theta < -theta_threshold) {
 805         *direction = UP;
 806     }
 807
 808     switch (*direction) {
 809     case RIGHT:
 810         *uf =  vec[2] / vec[0];
 811         *vf = -vec[1] / vec[0];
 812         break;
 813     case LEFT:
 814         *uf =  vec[2] / vec[0];
 815         *vf =  vec[1] / vec[0];
 816         break;
 817     case UP:
 818         *uf =  vec[0] / vec[1];
 819         *vf = -vec[2] / vec[1];
 820         break;
 821     case DOWN:
 822         *uf = -vec[0] / vec[1];
 823         *vf = -vec[2] / vec[1];
 824         break;
 825     case FRONT:
 826         *uf = -vec[0] / vec[2];
 827         *vf =  vec[1] / vec[2];
 828         break;
 829     case BACK:
 830         *uf = -vec[0] / vec[2];
 831         *vf = -vec[1] / vec[2];
 832         break;
 833     default:
 834         av_assert0(0);
 835     }
 836
 837     face = s->in_cubemap_face_order[*direction];
 838     rotate_cube_face(uf, vf, s->in_cubemap_face_rotation[face]);
 839
 840     (*uf) *= s->input_mirror_modifier[0];
 841     (*vf) *= s->input_mirror_modifier[1];
 842 }
 843
 844 /**
 845  * Find position on another cube face in case of overflow/underflow.
 846  * Used for calculation of interpolation window.
 847  *
 848  * @param s filter context
 849  * @param uf horizontal cubemap coordinate
 850  * @param vf vertical cubemap coordinate
 851  * @param direction direction of view
 852  * @param new_uf new horizontal cubemap coordinate
 853  * @param new_vf new vertical cubemap coordinate
 854  * @param face face position on cubemap
 855  */
 856 static void process_cube_coordinates(const V360Context *s,
 857                                      float uf, float vf, int direction,
 858                                      float *new_uf, float *new_vf, int *face)
 859 {
 860     /*
 861      *  Cubemap orientation
 862      *
 863      *           width
 864      *         <------->
 865      *         +-------+
 866      *         |       |                              U
 867      *         | up    |                   h       ------->
 868      * +-------+-------+-------+-------+ ^ e      |
 869      * |       |       |       |       | | i    V |
 870      * | left  | front | right | back  | | g      |
 871      * +-------+-------+-------+-------+ v h      v
 872      *         |       |                   t
 873      *         | down  |
 874      *         +-------+
 875      */
 876
 877     *face = s->in_cubemap_face_order[direction];
 878     rotate_cube_face_inverse(&uf, &vf, s->in_cubemap_face_rotation[*face]);
 879
 880     if ((uf < -1.f || uf >= 1.f) && (vf < -1.f || vf >= 1.f)) {
 881         // There are no pixels to use in this case
 882         *new_uf = uf;
 883         *new_vf = vf;
 884     } else if (uf < -1.f) {
 885         uf += 2.f;
 886         switch (direction) {
 887         case RIGHT:
 888             direction = FRONT;
 889             *new_uf =  uf;
 890             *new_vf =  vf;
 891             break;
 892         case LEFT:
 893             direction = BACK;
 894             *new_uf =  uf;
 895             *new_vf =  vf;
 896             break;
 897         case UP:
 898             direction = LEFT;
 899             *new_uf =  vf;
 900             *new_vf = -uf;
 901             break;
 902         case DOWN:
 903             direction = LEFT;
 904             *new_uf = -vf;
 905             *new_vf =  uf;
 906             break;
 907         case FRONT:
 908             direction = LEFT;
 909             *new_uf =  uf;
 910             *new_vf =  vf;
 911             break;
 912         case BACK:
 913             direction = RIGHT;
 914             *new_uf =  uf;
 915             *new_vf =  vf;
 916             break;
 917         default:
 918             av_assert0(0);
 919         }
 920     } else if (uf >= 1.f) {
 921         uf -= 2.f;
 922         switch (direction) {
 923         case RIGHT:
 924             direction = BACK;
 925             *new_uf =  uf;
 926             *new_vf =  vf;
 927             break;
 928         case LEFT:
 929             direction = FRONT;
 930             *new_uf =  uf;
 931             *new_vf =  vf;
 932             break;
 933         case UP:
 934             direction = RIGHT;
 935             *new_uf = -vf;
 936             *new_vf =  uf;
 937             break;
 938         case DOWN:
 939             direction = RIGHT;
 940             *new_uf =  vf;
 941             *new_vf = -uf;
 942             break;
 943         case FRONT:
 944             direction = RIGHT;
 945             *new_uf =  uf;
 946             *new_vf =  vf;
 947             break;
 948         case BACK:
 949             direction = LEFT;
 950             *new_uf =  uf;
 951             *new_vf =  vf;
 952             break;
 953         default:
 954             av_assert0(0);
 955         }
 956     } else if (vf < -1.f) {
 957         vf += 2.f;
 958         switch (direction) {
 959         case RIGHT:
 960             direction = UP;
 961             *new_uf =  vf;
 962             *new_vf = -uf;
 963             break;
 964         case LEFT:
 965             direction = UP;
 966             *new_uf = -vf;
 967             *new_vf =  uf;
 968             break;
 969         case UP:
 970             direction = BACK;
 971             *new_uf = -uf;
 972             *new_vf = -vf;
 973             break;
 974         case DOWN:
 975             direction = FRONT;
 976             *new_uf =  uf;
 977             *new_vf =  vf;
 978             break;
 979         case FRONT:
 980             direction = UP;
 981             *new_uf =  uf;
 982             *new_vf =  vf;
 983             break;
 984         case BACK:
 985             direction = UP;
 986             *new_uf = -uf;
 987             *new_vf = -vf;
 988             break;
 989         default:
 990             av_assert0(0);
 991         }
 992     } else if (vf >= 1.f) {
 993         vf -= 2.f;
 994         switch (direction) {
 995         case RIGHT:
 996             direction = DOWN;
 997             *new_uf = -vf;
 998             *new_vf =  uf;
 999             break;
1000         case LEFT:
1001             direction = DOWN;
1002             *new_uf =  vf;
1003             *new_vf = -uf;
1004             break;
1005         case UP:
1006             direction = FRONT;
1007             *new_uf =  uf;
1008             *new_vf =  vf;
1009             break;
1010         case DOWN:
1011             direction = BACK;
1012             *new_uf = -uf;
1013             *new_vf = -vf;
1014             break;
1015         case FRONT:
1016             direction = DOWN;
1017             *new_uf =  uf;
1018             *new_vf =  vf;
1019             break;
1020         case BACK:
1021             direction = DOWN;
1022             *new_uf = -uf;
1023             *new_vf = -vf;
1024             break;
1025         default:
1026             av_assert0(0);
1027         }
1028     } else {
1029         // Inside cube face
1030         *new_uf = uf;
1031         *new_vf = vf;
1032     }
1033
1034     *face = s->in_cubemap_face_order[direction];
1035     rotate_cube_face(new_uf, new_vf, s->in_cubemap_face_rotation[*face]);
1036 }
1037
1038 /**
1039  * Calculate 3D coordinates on sphere for corresponding frame position in cubemap3x2 format.
1040  *
1041  * @param s filter context
1042  * @param i horizontal position on frame [0, height)
1043  * @param j vertical position on frame [0, width)
1044  * @param width frame width
1045  * @param height frame height
1046  * @param vec coordinates on sphere
1047  */
1048 static void cube3x2_to_xyz(const V360Context *s,
1049                            int i, int j, int width, int height,
1050                            float *vec)
1051 {
1052     const float ew = width  / 3.f;
1053     const float eh = height / 2.f;
1054
1055     const int u_face = floorf(i / ew);
1056     const int v_face = floorf(j / eh);
1057     const int face = u_face + 3 * v_face;
1058
1059     const int u_shift = ceilf(ew * u_face);
1060     const int v_shift = ceilf(eh * v_face);
1061     const int ewi = ceilf(ew * (u_face + 1)) - u_shift;
1062     const int ehi = ceilf(eh * (v_face + 1)) - v_shift;
1063
1064     const float uf = 2.f * (i - u_shift) / ewi - 1.f;
1065     const float vf = 2.f * (j - v_shift) / ehi - 1.f;
1066
1067     cube_to_xyz(s, uf, vf, face, vec);
1068 }
1069
1070 /**
1071  * Calculate frame position in cubemap3x2 format for corresponding 3D coordinates on sphere.
1072  *
1073  * @param s filter context
1074  * @param vec coordinates on sphere
1075  * @param width frame width
1076  * @param height frame height
1077  * @param us horizontal coordinates for interpolation window
1078  * @param vs vertical coordinates for interpolation window
1079  * @param du horizontal relative coordinate
1080  * @param dv vertical relative coordinate
1081  */
1082 static void xyz_to_cube3x2(const V360Context *s,
1083                            const float *vec, int width, int height,
1084                            uint16_t us[4][4], uint16_t vs[4][4], float *du, float *dv)
1085 {
1086     const float ew = width  / 3.f;
1087     const float eh = height / 2.f;
1088     float uf, vf;
1089     int ui, vi;
1090     int ewi, ehi;
1091     int i, j;
1092     int direction, face;
1093     int u_face, v_face;
1094
1095     xyz_to_cube(s, vec, &uf, &vf, &direction);
1096
1097     uf *= (1.f - s->in_pad);
1098     vf *= (1.f - s->in_pad);
1099
1100     face = s->in_cubemap_face_order[direction];
1101     u_face = face % 3;
1102     v_face = face / 3;
1103     ewi = ceilf(ew * (u_face + 1)) - ceilf(ew * u_face);
1104     ehi = ceilf(eh * (v_face + 1)) - ceilf(eh * v_face);
1105
1106     uf = 0.5f * ewi * (uf + 1.f);
1107     vf = 0.5f * ehi * (vf + 1.f);
1108
1109     ui = floorf(uf);
1110     vi = floorf(vf);
1111
1112     *du = uf - ui;
1113     *dv = vf - vi;
1114
1115     for (i = -1; i < 3; i++) {
1116         for (j = -1; j < 3; j++) {
1117             int new_ui = ui + j;
1118             int new_vi = vi + i;
1119             int u_shift, v_shift;
1120             int new_ewi, new_ehi;
1121
1122             if (new_ui >= 0 && new_ui < ewi && new_vi >= 0 && new_vi < ehi) {
1123                 face = s->in_cubemap_face_order[direction];
1124
1125                 u_face = face % 3;
1126                 v_face = face / 3;
1127                 u_shift = ceilf(ew * u_face);
1128                 v_shift = ceilf(eh * v_face);
1129             } else {
1130                 uf = 2.f * new_ui / ewi - 1.f;
1131                 vf = 2.f * new_vi / ehi - 1.f;
1132
1133                 uf /= (1.f - s->in_pad);
1134                 vf /= (1.f - s->in_pad);
1135
1136                 process_cube_coordinates(s, uf, vf, direction, &uf, &vf, &face);
1137
1138                 uf *= (1.f - s->in_pad);
1139                 vf *= (1.f - s->in_pad);
1140
1141                 u_face = face % 3;
1142                 v_face = face / 3;
1143                 u_shift = ceilf(ew * u_face);
1144                 v_shift = ceilf(eh * v_face);
1145                 new_ewi = ceilf(ew * (u_face + 1)) - u_shift;
1146                 new_ehi = ceilf(eh * (v_face + 1)) - v_shift;
1147
1148                 new_ui = av_clip(roundf(0.5f * new_ewi * (uf + 1.f)), 0, new_ewi - 1);
1149                 new_vi = av_clip(roundf(0.5f * new_ehi * (vf + 1.f)), 0, new_ehi - 1);
1150             }
1151
1152             us[i + 1][j + 1] = u_shift + new_ui;
1153             vs[i + 1][j + 1] = v_shift + new_vi;
1154         }
1155     }
1156 }
1157
1158 /**
1159  * Calculate 3D coordinates on sphere for corresponding frame position in cubemap1x6 format.
1160  *
1161  * @param s filter context
1162  * @param i horizontal position on frame [0, height)
1163  * @param j vertical position on frame [0, width)
1164  * @param width frame width
1165  * @param height frame height
1166  * @param vec coordinates on sphere
1167  */
1168 static void cube1x6_to_xyz(const V360Context *s,
1169                            int i, int j, int width, int height,
1170                            float *vec)
1171 {
1172     const float ew = width;
1173     const float eh = height / 6.f;
1174
1175     const int face = floorf(j / eh);
1176
1177     const int v_shift = ceilf(eh * face);
1178     const int ehi = ceilf(eh * (face + 1)) - v_shift;
1179
1180     const float uf = 2.f *  i            / ew  - 1.f;
1181     const float vf = 2.f * (j - v_shift) / ehi - 1.f;
1182
1183     cube_to_xyz(s, uf, vf, face, vec);
1184 }
1185
1186 /**
1187  * Calculate 3D coordinates on sphere for corresponding frame position in cubemap6x1 format.
1188  *
1189  * @param s filter context
1190  * @param i horizontal position on frame [0, height)
1191  * @param j vertical position on frame [0, width)
1192  * @param width frame width
1193  * @param height frame height
1194  * @param vec coordinates on sphere
1195  */
1196 static void cube6x1_to_xyz(const V360Context *s,
1197                            int i, int j, int width, int height,
1198                            float *vec)
1199 {
1200     const float ew = width / 6.f;
1201     const float eh = height;
1202
1203     const int face = floorf(i / ew);
1204
1205     const int u_shift = ceilf(ew * face);
1206     const int ewi = ceilf(ew * (face + 1)) - u_shift;
1207
1208     const float uf = 2.f * (i - u_shift) / ewi - 1.f;
1209     const float vf = 2.f *  j            / eh  - 1.f;
1210
1211     cube_to_xyz(s, uf, vf, face, vec);
1212 }
1213
1214 /**
1215  * Calculate frame position in cubemap1x6 format for corresponding 3D coordinates on sphere.
1216  *
1217  * @param s filter context
1218  * @param vec coordinates on sphere
1219  * @param width frame width
1220  * @param height frame height
1221  * @param us horizontal coordinates for interpolation window
1222  * @param vs vertical coordinates for interpolation window
1223  * @param du horizontal relative coordinate
1224  * @param dv vertical relative coordinate
1225  */
1226 static void xyz_to_cube1x6(const V360Context *s,
1227                            const float *vec, int width, int height,
1228                            uint16_t us[4][4], uint16_t vs[4][4], float *du, float *dv)
1229 {
1230     const float eh = height / 6.f;
1231     const int ewi = width;
1232     float uf, vf;
1233     int ui, vi;
1234     int ehi;
1235     int i, j;
1236     int direction, face;
1237
1238     xyz_to_cube(s, vec, &uf, &vf, &direction);
1239
1240     uf *= (1.f - s->in_pad);
1241     vf *= (1.f - s->in_pad);
1242
1243     face = s->in_cubemap_face_order[direction];
1244     ehi = ceilf(eh * (face + 1)) - ceilf(eh * face);
1245
1246     uf = 0.5f * ewi * (uf + 1.f);
1247     vf = 0.5f * ehi * (vf + 1.f);
1248
1249     ui = floorf(uf);
1250     vi = floorf(vf);
1251
1252     *du = uf - ui;
1253     *dv = vf - vi;
1254
1255     for (i = -1; i < 3; i++) {
1256         for (j = -1; j < 3; j++) {
1257             int new_ui = ui + j;
1258             int new_vi = vi + i;
1259             int v_shift;
1260             int new_ehi;
1261
1262             if (new_ui >= 0 && new_ui < ewi && new_vi >= 0 && new_vi < ehi) {
1263                 face = s->in_cubemap_face_order[direction];
1264
1265                 v_shift = ceilf(eh * face);
1266             } else {
1267                 uf = 2.f * new_ui / ewi - 1.f;
1268                 vf = 2.f * new_vi / ehi - 1.f;
1269
1270                 uf /= (1.f - s->in_pad);
1271                 vf /= (1.f - s->in_pad);
1272
1273                 process_cube_coordinates(s, uf, vf, direction, &uf, &vf, &face);
1274
1275                 uf *= (1.f - s->in_pad);
1276                 vf *= (1.f - s->in_pad);
1277
1278                 v_shift = ceilf(eh * face);
1279                 new_ehi = ceilf(eh * (face + 1)) - v_shift;
1280
1281                 new_ui = av_clip(roundf(0.5f *     ewi * (uf + 1.f)), 0,     ewi - 1);
1282                 new_vi = av_clip(roundf(0.5f * new_ehi * (vf + 1.f)), 0, new_ehi - 1);
1283             }
1284
1285             us[i + 1][j + 1] =           new_ui;
1286             vs[i + 1][j + 1] = v_shift + new_vi;
1287         }
1288     }
1289 }
1290
1291 /**
1292  * Calculate frame position in cubemap6x1 format for corresponding 3D coordinates on sphere.
1293  *
1294  * @param s filter context
1295  * @param vec coordinates on sphere
1296  * @param width frame width
1297  * @param height frame height
1298  * @param us horizontal coordinates for interpolation window
1299  * @param vs vertical coordinates for interpolation window
1300  * @param du horizontal relative coordinate
1301  * @param dv vertical relative coordinate
1302  */
1303 static void xyz_to_cube6x1(const V360Context *s,
1304                            const float *vec, int width, int height,
1305                            uint16_t us[4][4], uint16_t vs[4][4], float *du, float *dv)
1306 {
1307     const float ew = width / 6.f;
1308     const int ehi = height;
1309     float uf, vf;
1310     int ui, vi;
1311     int ewi;
1312     int i, j;
1313     int direction, face;
1314
1315     xyz_to_cube(s, vec, &uf, &vf, &direction);
1316
1317     uf *= (1.f - s->in_pad);
1318     vf *= (1.f - s->in_pad);
1319
1320     face = s->in_cubemap_face_order[direction];
1321     ewi = ceilf(ew * (face + 1)) - ceilf(ew * face);
1322
1323     uf = 0.5f * ewi * (uf + 1.f);
1324     vf = 0.5f * ehi * (vf + 1.f);
1325
1326     ui = floorf(uf);
1327     vi = floorf(vf);
1328
1329     *du = uf - ui;
1330     *dv = vf - vi;
1331
1332     for (i = -1; i < 3; i++) {
1333         for (j = -1; j < 3; j++) {
1334             int new_ui = ui + j;
1335             int new_vi = vi + i;
1336             int u_shift;
1337             int new_ewi;
1338
1339             if (new_ui >= 0 && new_ui < ewi && new_vi >= 0 && new_vi < ehi) {
1340                 face = s->in_cubemap_face_order[direction];
1341
1342                 u_shift = ceilf(ew * face);
1343             } else {
1344                 uf = 2.f * new_ui / ewi - 1.f;
1345                 vf = 2.f * new_vi / ehi - 1.f;
1346
1347                 uf /= (1.f - s->in_pad);
1348                 vf /= (1.f - s->in_pad);
1349
1350                 process_cube_coordinates(s, uf, vf, direction, &uf, &vf, &face);
1351
1352                 uf *= (1.f - s->in_pad);
1353                 vf *= (1.f - s->in_pad);
1354
1355                 u_shift = ceilf(ew * face);
1356                 new_ewi = ceilf(ew * (face + 1)) - u_shift;
1357
1358                 new_ui = av_clip(roundf(0.5f * new_ewi * (uf + 1.f)), 0, new_ewi - 1);
1359                 new_vi = av_clip(roundf(0.5f *     ehi * (vf + 1.f)), 0,     ehi - 1);
1360             }
1361
1362             us[i + 1][j + 1] = u_shift + new_ui;
1363             vs[i + 1][j + 1] =           new_vi;
1364         }
1365     }
1366 }
1367
1368 /**
1369  * Calculate 3D coordinates on sphere for corresponding frame position in equirectangular format.
1370  *
1371  * @param s filter context
1372  * @param i horizontal position on frame [0, height)
1373  * @param j vertical position on frame [0, width)
1374  * @param width frame width
1375  * @param height frame height
1376  * @param vec coordinates on sphere
1377  */
1378 static void equirect_to_xyz(const V360Context *s,
1379                             int i, int j, int width, int height,
1380                             float *vec)
1381 {
1382     const float phi   = ((2.f * i) / width  - 1.f) * M_PI;
1383     const float theta = ((2.f * j) / height - 1.f) * M_PI_2;
1384
1385     const float sin_phi   = sinf(phi);
1386     const float cos_phi   = cosf(phi);
1387     const float sin_theta = sinf(theta);
1388     const float cos_theta = cosf(theta);
1389
1390     vec[0] =  cos_theta * sin_phi;
1391     vec[1] = -sin_theta;
1392     vec[2] = -cos_theta * cos_phi;
1393 }
1394
1395 /**
1396  * Calculate 3D coordinates on sphere for corresponding frame position in stereographic format.
1397  *
1398  * @param s filter context
1399  * @param i horizontal position on frame [0, height)
1400  * @param j vertical position on frame [0, width)
1401  * @param width frame width
1402  * @param height frame height
1403  * @param vec coordinates on sphere
1404  */
1405 static void stereographic_to_xyz(const V360Context *s,
1406                                  int i, int j, int width, int height,
1407                                  float *vec)
1408 {
1409     const float x = ((2.f * i) / width  - 1.f) * (s->h_fov / 180.f) * M_PI;
1410     const float y = ((2.f * j) / height - 1.f) * (s->v_fov /  90.f) * M_PI_2;
1411     const float xy = x * x + y * y;
1412
1413     vec[0] = 2.f * x / (1.f + xy);
1414     vec[1] = (-1.f + xy) / (1.f + xy);
1415     vec[2] = 2.f * y / (1.f + xy);
1416
1417     normalize_vector(vec);
1418 }
1419
1420 /**
1421  * Calculate frame position in equirectangular format for corresponding 3D coordinates on sphere.
1422  *
1423  * @param s filter context
1424  * @param vec coordinates on sphere
1425  * @param width frame width
1426  * @param height frame height
1427  * @param us horizontal coordinates for interpolation window
1428  * @param vs vertical coordinates for interpolation window
1429  * @param du horizontal relative coordinate
1430  * @param dv vertical relative coordinate
1431  */
1432 static void xyz_to_equirect(const V360Context *s,
1433                             const float *vec, int width, int height,
1434                             uint16_t us[4][4], uint16_t vs[4][4], float *du, float *dv)
1435 {
1436     const float phi   = atan2f(vec[0], -vec[2]) * s->input_mirror_modifier[0];
1437     const float theta = asinf(-vec[1]) * s->input_mirror_modifier[1];
1438     float uf, vf;
1439     int ui, vi;
1440     int i, j;
1441
1442     uf = (phi   / M_PI   + 1.f) * width  / 2.f;
1443     vf = (theta / M_PI_2 + 1.f) * height / 2.f;
1444     ui = floorf(uf);
1445     vi = floorf(vf);
1446
1447     *du = uf - ui;
1448     *dv = vf - vi;
1449
1450     for (i = -1; i < 3; i++) {
1451         for (j = -1; j < 3; j++) {
1452             us[i + 1][j + 1] = mod(ui + j, width);
1453             vs[i + 1][j + 1] = av_clip(vi + i, 0, height - 1);
1454         }
1455     }
1456 }
1457
1458 /**
1459  * Prepare data for processing equi-angular cubemap input format.
1460  *
1461  * @param ctx filter context
1462
1463  * @return error code
1464  */
1465 static int prepare_eac_in(AVFilterContext *ctx)
1466 {
1467     V360Context *s = ctx->priv;
1468
1469     if (s->ih_flip && s->iv_flip) {
1470         s->in_cubemap_face_order[RIGHT] = BOTTOM_LEFT;
1471         s->in_cubemap_face_order[LEFT]  = BOTTOM_RIGHT;
1472         s->in_cubemap_face_order[UP]    = TOP_LEFT;
1473         s->in_cubemap_face_order[DOWN]  = TOP_RIGHT;
1474         s->in_cubemap_face_order[FRONT] = BOTTOM_MIDDLE;
1475         s->in_cubemap_face_order[BACK]  = TOP_MIDDLE;
1476     } else if (s->ih_flip) {
1477         s->in_cubemap_face_order[RIGHT] = TOP_LEFT;
1478         s->in_cubemap_face_order[LEFT]  = TOP_RIGHT;
1479         s->in_cubemap_face_order[UP]    = BOTTOM_LEFT;
1480         s->in_cubemap_face_order[DOWN]  = BOTTOM_RIGHT;
1481         s->in_cubemap_face_order[FRONT] = TOP_MIDDLE;
1482         s->in_cubemap_face_order[BACK]  = BOTTOM_MIDDLE;
1483     } else if (s->iv_flip) {
1484         s->in_cubemap_face_order[RIGHT] = BOTTOM_RIGHT;
1485         s->in_cubemap_face_order[LEFT]  = BOTTOM_LEFT;
1486         s->in_cubemap_face_order[UP]    = TOP_RIGHT;
1487         s->in_cubemap_face_order[DOWN]  = TOP_LEFT;
1488         s->in_cubemap_face_order[FRONT] = BOTTOM_MIDDLE;
1489         s->in_cubemap_face_order[BACK]  = TOP_MIDDLE;
1490     } else {
1491         s->in_cubemap_face_order[RIGHT] = TOP_RIGHT;
1492         s->in_cubemap_face_order[LEFT]  = TOP_LEFT;
1493         s->in_cubemap_face_order[UP]    = BOTTOM_RIGHT;
1494         s->in_cubemap_face_order[DOWN]  = BOTTOM_LEFT;
1495         s->in_cubemap_face_order[FRONT] = TOP_MIDDLE;
1496         s->in_cubemap_face_order[BACK]  = BOTTOM_MIDDLE;
1497     }
1498
1499     if (s->iv_flip) {
1500         s->in_cubemap_face_rotation[TOP_LEFT]      = ROT_270;
1501         s->in_cubemap_face_rotation[TOP_MIDDLE]    = ROT_90;
1502         s->in_cubemap_face_rotation[TOP_RIGHT]     = ROT_270;
1503         s->in_cubemap_face_rotation[BOTTOM_LEFT]   = ROT_0;
1504         s->in_cubemap_face_rotation[BOTTOM_MIDDLE] = ROT_0;
1505         s->in_cubemap_face_rotation[BOTTOM_RIGHT]  = ROT_0;
1506     } else {
1507         s->in_cubemap_face_rotation[TOP_LEFT]      = ROT_0;
1508         s->in_cubemap_face_rotation[TOP_MIDDLE]    = ROT_0;
1509         s->in_cubemap_face_rotation[TOP_RIGHT]     = ROT_0;
1510         s->in_cubemap_face_rotation[BOTTOM_LEFT]   = ROT_270;
1511         s->in_cubemap_face_rotation[BOTTOM_MIDDLE] = ROT_90;
1512         s->in_cubemap_face_rotation[BOTTOM_RIGHT]  = ROT_270;
1513     }
1514
1515     return 0;
1516 }
1517
1518 /**
1519  * Prepare data for processing equi-angular cubemap output format.
1520  *
1521  * @param ctx filter context
1522  *
1523  * @return error code
1524  */
1525 static int prepare_eac_out(AVFilterContext *ctx)
1526 {
1527     V360Context *s = ctx->priv;
1528
1529     s->out_cubemap_direction_order[TOP_LEFT]      = LEFT;
1530     s->out_cubemap_direction_order[TOP_MIDDLE]    = FRONT;
1531     s->out_cubemap_direction_order[TOP_RIGHT]     = RIGHT;
1532     s->out_cubemap_direction_order[BOTTOM_LEFT]   = DOWN;
1533     s->out_cubemap_direction_order[BOTTOM_MIDDLE] = BACK;
1534     s->out_cubemap_direction_order[BOTTOM_RIGHT]  = UP;
1535
1536     s->out_cubemap_face_rotation[TOP_LEFT]      = ROT_0;
1537     s->out_cubemap_face_rotation[TOP_MIDDLE]    = ROT_0;
1538     s->out_cubemap_face_rotation[TOP_RIGHT]     = ROT_0;
1539     s->out_cubemap_face_rotation[BOTTOM_LEFT]   = ROT_270;
1540     s->out_cubemap_face_rotation[BOTTOM_MIDDLE] = ROT_90;
1541     s->out_cubemap_face_rotation[BOTTOM_RIGHT]  = ROT_270;
1542
1543     return 0;
1544 }
1545
1546 /**
1547  * Calculate 3D coordinates on sphere for corresponding frame position in equi-angular cubemap format.
1548  *
1549  * @param s filter context
1550  * @param i horizontal position on frame [0, height)
1551  * @param j vertical position on frame [0, width)
1552  * @param width frame width
1553  * @param height frame height
1554  * @param vec coordinates on sphere
1555  */
1556 static void eac_to_xyz(const V360Context *s,
1557                        int i, int j, int width, int height,
1558                        float *vec)
1559 {
1560     const float pixel_pad = 2;
1561     const float u_pad = pixel_pad / width;
1562     const float v_pad = pixel_pad / height;
1563
1564     int u_face, v_face, face;
1565
1566     float l_x, l_y, l_z;
1567
1568     float uf = (float)i / width;
1569     float vf = (float)j / height;
1570
1571     // EAC has 2-pixel padding on faces except between faces on the same row
1572     // Padding pixels seems not to be stretched with tangent as regular pixels
1573     // Formulas below approximate original padding as close as I could get experimentally
1574
1575     // Horizontal padding
1576     uf = 3.f * (uf - u_pad) / (1.f - 2.f * u_pad);
1577     if (uf < 0.f) {
1578         u_face = 0;
1579         uf -= 0.5f;
1580     } else if (uf >= 3.f) {
1581         u_face = 2;
1582         uf -= 2.5f;
1583     } else {
1584         u_face = floorf(uf);
1585         uf = fmodf(uf, 1.f) - 0.5f;
1586     }
1587
1588     // Vertical padding
1589     v_face = floorf(vf * 2.f);
1590     vf = (vf - v_pad - 0.5f * v_face) / (0.5f - 2.f * v_pad) - 0.5f;
1591
1592     if (uf >= -0.5f && uf < 0.5f) {
1593         uf = tanf(M_PI_2 * uf);
1594     } else {
1595         uf = 2.f * uf;
1596     }
1597     if (vf >= -0.5f && vf < 0.5f) {
1598         vf = tanf(M_PI_2 * vf);
1599     } else {
1600         vf = 2.f * vf;
1601     }
1602
1603     face = u_face + 3 * v_face;
1604
1605     switch (face) {
1606     case TOP_LEFT:
1607         l_x = -1.f;
1608         l_y = -vf;
1609         l_z = -uf;
1610         break;
1611     case TOP_MIDDLE:
1612         l_x =  uf;
1613         l_y = -vf;
1614         l_z = -1.f;
1615         break;
1616     case TOP_RIGHT:
1617         l_x =  1.f;
1618         l_y = -vf;
1619         l_z =  uf;
1620         break;
1621     case BOTTOM_LEFT:
1622         l_x = -vf;
1623         l_y = -1.f;
1624         l_z =  uf;
1625         break;
1626     case BOTTOM_MIDDLE:
1627         l_x = -vf;
1628         l_y =  uf;
1629         l_z =  1.f;
1630         break;
1631     case BOTTOM_RIGHT:
1632         l_x = -vf;
1633         l_y =  1.f;
1634         l_z = -uf;
1635         break;
1636     default:
1637         av_assert0(0);
1638     }
1639
1640     vec[0] = l_x;
1641     vec[1] = l_y;
1642     vec[2] = l_z;
1643
1644     normalize_vector(vec);
1645 }
1646
1647 /**
1648  * Calculate frame position in equi-angular cubemap format for corresponding 3D coordinates on sphere.
1649  *
1650  * @param s filter context
1651  * @param vec coordinates on sphere
1652  * @param width frame width
1653  * @param height frame height
1654  * @param us horizontal coordinates for interpolation window
1655  * @param vs vertical coordinates for interpolation window
1656  * @param du horizontal relative coordinate
1657  * @param dv vertical relative coordinate
1658  */
1659 static void xyz_to_eac(const V360Context *s,
1660                        const float *vec, int width, int height,
1661                        uint16_t us[4][4], uint16_t vs[4][4], float *du, float *dv)
1662 {
1663     const float pixel_pad = 2;
1664     const float u_pad = pixel_pad / width;
1665     const float v_pad = pixel_pad / height;
1666
1667     float uf, vf;
1668     int ui, vi;
1669     int i, j;
1670     int direction, face;
1671     int u_face, v_face;
1672
1673     xyz_to_cube(s, vec, &uf, &vf, &direction);
1674
1675     face = s->in_cubemap_face_order[direction];
1676     u_face = face % 3;
1677     v_face = face / 3;
1678
1679     uf = M_2_PI * atanf(uf) + 0.5f;
1680     vf = M_2_PI * atanf(vf) + 0.5f;
1681
1682     // These formulas are inversed from eac_to_xyz ones
1683     uf = (uf + u_face) * (1.f - 2.f * u_pad) / 3.f + u_pad;
1684     vf = vf * (0.5f - 2.f * v_pad) + v_pad + 0.5f * v_face;
1685
1686     uf *= width;
1687     vf *= height;
1688
1689     ui = floorf(uf);
1690     vi = floorf(vf);
1691
1692     *du = uf - ui;
1693     *dv = vf - vi;
1694
1695     for (i = -1; i < 3; i++) {
1696         for (j = -1; j < 3; j++) {
1697             us[i + 1][j + 1] = av_clip(ui + j, 0, width  - 1);
1698             vs[i + 1][j + 1] = av_clip(vi + i, 0, height - 1);
1699         }
1700     }
1701 }
1702
1703 /**
1704  * Prepare data for processing flat output format.
1705  *
1706  * @param ctx filter context
1707  *
1708  * @return error code
1709  */
1710 static int prepare_flat_out(AVFilterContext *ctx)
1711 {
1712     V360Context *s = ctx->priv;
1713
1714     const float h_angle = 0.5f * s->h_fov * M_PI / 180.f;
1715     const float v_angle = 0.5f * s->v_fov * M_PI / 180.f;
1716
1717     const float sin_phi   = sinf(h_angle);
1718     const float cos_phi   = cosf(h_angle);
1719     const float sin_theta = sinf(v_angle);
1720     const float cos_theta = cosf(v_angle);
1721
1722     s->flat_range[0] =  cos_theta * sin_phi;
1723     s->flat_range[1] =  sin_theta;
1724     s->flat_range[2] = -cos_theta * cos_phi;
1725
1726     return 0;
1727 }
1728
1729 /**
1730  * Calculate 3D coordinates on sphere for corresponding frame position in flat format.
1731  *
1732  * @param s filter context
1733  * @param i horizontal position on frame [0, height)
1734  * @param j vertical position on frame [0, width)
1735  * @param width frame width
1736  * @param height frame height
1737  * @param vec coordinates on sphere
1738  */
1739 static void flat_to_xyz(const V360Context *s,
1740                         int i, int j, int width, int height,
1741                         float *vec)
1742 {
1743     const float l_x =  s->flat_range[0] * (2.f * i / width  - 1.f);
1744     const float l_y = -s->flat_range[1] * (2.f * j / height - 1.f);
1745     const float l_z =  s->flat_range[2];
1746
1747     vec[0] = l_x;
1748     vec[1] = l_y;
1749     vec[2] = l_z;
1750
1751     normalize_vector(vec);
1752 }
1753
1754 /**
1755  * Calculate frame position in dual fisheye format for corresponding 3D coordinates on sphere.
1756  *
1757  * @param s filter context
1758  * @param vec coordinates on sphere
1759  * @param width frame width
1760  * @param height frame height
1761  * @param us horizontal coordinates for interpolation window
1762  * @param vs vertical coordinates for interpolation window
1763  * @param du horizontal relative coordinate
1764  * @param dv vertical relative coordinate
1765  */
1766 static void xyz_to_dfisheye(const V360Context *s,
1767                             const float *vec, int width, int height,
1768                             uint16_t us[4][4], uint16_t vs[4][4], float *du, float *dv)
1769 {
1770     const float scale = 1.f - s->in_pad;
1771
1772     const float ew = width / 2.f;
1773     const float eh = height;
1774
1775     const float phi   = atan2f(-vec[1], -vec[0]) * s->input_mirror_modifier[0];
1776     const float theta = acosf(fabsf(vec[2])) / M_PI * s->input_mirror_modifier[1];
1777
1778     float uf = (theta * cosf(phi) * scale + 0.5f) * ew;
1779     float vf = (theta * sinf(phi) * scale + 0.5f) * eh;
1780
1781     int ui, vi;
1782     int u_shift;
1783     int i, j;
1784
1785     if (vec[2] >= 0) {
1786         u_shift = 0;
1787     } else {
1788         u_shift = ceilf(ew);
1789         uf = ew - uf;
1790     }
1791
1792     ui = floorf(uf);
1793     vi = floorf(vf);
1794
1795     *du = uf - ui;
1796     *dv = vf - vi;
1797
1798     for (i = -1; i < 3; i++) {
1799         for (j = -1; j < 3; j++) {
1800             us[i + 1][j + 1] = av_clip(u_shift + ui + j, 0, width  - 1);
1801             vs[i + 1][j + 1] = av_clip(          vi + i, 0, height - 1);
1802         }
1803     }
1804 }
1805
1806 /**
1807  * Calculate 3D coordinates on sphere for corresponding frame position in barrel facebook's format.
1808  *
1809  * @param s filter context
1810  * @param i horizontal position on frame [0, height)
1811  * @param j vertical position on frame [0, width)
1812  * @param width frame width
1813  * @param height frame height
1814  * @param vec coordinates on sphere
1815  */
1816 static void barrel_to_xyz(const V360Context *s,
1817                           int i, int j, int width, int height,
1818                           float *vec)
1819 {
1820     const float scale = 0.99f;
1821     float l_x, l_y, l_z;
1822
1823     if (i < 4 * width / 5) {
1824         const float theta_range = M_PI / 4.f;
1825
1826         const int ew = 4 * width / 5;
1827         const int eh = height;
1828
1829         const float phi   = ((2.f * i) / ew - 1.f) * M_PI        / scale;
1830         const float theta = ((2.f * j) / eh - 1.f) * theta_range / scale;
1831
1832         const float sin_phi   = sinf(phi);
1833         const float cos_phi   = cosf(phi);
1834         const float sin_theta = sinf(theta);
1835         const float cos_theta = cosf(theta);
1836
1837         l_x =  cos_theta * sin_phi;
1838         l_y = -sin_theta;
1839         l_z = -cos_theta * cos_phi;
1840     } else {
1841         const int ew = width  / 5;
1842         const int eh = height / 2;
1843
1844         float uf, vf;
1845
1846         if (j < eh) {   // UP
1847             uf = 2.f * (i - 4 * ew) / ew  - 1.f;
1848             vf = 2.f * (j         ) / eh - 1.f;
1849
1850             uf /= scale;
1851             vf /= scale;
1852
1853             l_x =  uf;
1854             l_y =  1.f;
1855             l_z = -vf;
1856         } else {            // DOWN
1857             uf = 2.f * (i - 4 * ew) / ew - 1.f;
1858             vf = 2.f * (j -     eh) / eh - 1.f;
1859
1860             uf /= scale;
1861             vf /= scale;
1862
1863             l_x =  uf;
1864             l_y = -1.f;
1865             l_z =  vf;
1866         }
1867     }
1868
1869     vec[0] = l_x;
1870     vec[1] = l_y;
1871     vec[2] = l_z;
1872
1873     normalize_vector(vec);
1874 }
1875
1876 /**
1877  * Calculate frame position in barrel facebook's format for corresponding 3D coordinates on sphere.
1878  *
1879  * @param s filter context
1880  * @param vec coordinates on sphere
1881  * @param width frame width
1882  * @param height frame height
1883  * @param us horizontal coordinates for interpolation window
1884  * @param vs vertical coordinates for interpolation window
1885  * @param du horizontal relative coordinate
1886  * @param dv vertical relative coordinate
1887  */
1888 static void xyz_to_barrel(const V360Context *s,
1889                           const float *vec, int width, int height,
1890                           uint16_t us[4][4], uint16_t vs[4][4], float *du, float *dv)
1891 {
1892     const float scale = 0.99f;
1893
1894     const float phi   = atan2f(vec[0], -vec[2]) * s->input_mirror_modifier[0];
1895     const float theta = asinf(-vec[1]) * s->input_mirror_modifier[1];
1896     const float theta_range = M_PI / 4.f;
1897
1898     int ew, eh;
1899     int u_shift, v_shift;
1900     float uf, vf;
1901     int ui, vi;
1902     int i, j;
1903
1904     if (theta > -theta_range && theta < theta_range) {
1905         ew = 4 * width / 5;
1906         eh = height;
1907
1908         u_shift = s->ih_flip ? width / 5 : 0;
1909         v_shift = 0;
1910
1911         uf = (phi   / M_PI        * scale + 1.f) * ew / 2.f;
1912         vf = (theta / theta_range * scale + 1.f) * eh / 2.f;
1913     } else {
1914         ew = width  / 5;
1915         eh = height / 2;
1916
1917         u_shift = s->ih_flip ? 0 : 4 * ew;
1918
1919         if (theta < 0.f) {  // UP
1920             uf =  vec[0] / vec[1];
1921             vf = -vec[2] / vec[1];
1922             v_shift = 0;
1923         } else {            // DOWN
1924             uf = -vec[0] / vec[1];
1925             vf = -vec[2] / vec[1];
1926             v_shift = eh;
1927         }
1928
1929         uf *= s->input_mirror_modifier[0] * s->input_mirror_modifier[1];
1930         vf *= s->input_mirror_modifier[1];
1931
1932         uf = 0.5f * ew * (uf * scale + 1.f);
1933         vf = 0.5f * eh * (vf * scale + 1.f);
1934     }
1935
1936     ui = floorf(uf);
1937     vi = floorf(vf);
1938
1939     *du = uf - ui;
1940     *dv = vf - vi;
1941
1942     for (i = -1; i < 3; i++) {
1943         for (j = -1; j < 3; j++) {
1944             us[i + 1][j + 1] = u_shift + av_clip(ui + j, 0, ew - 1);
1945             vs[i + 1][j + 1] = v_shift + av_clip(vi + i, 0, eh - 1);
1946         }
1947     }
1948 }
1949
1950 static void multiply_matrix(float c[3][3], const float a[3][3], const float b[3][3])
1951 {
1952     for (int i = 0; i < 3; i++) {
1953         for (int j = 0; j < 3; j++) {
1954             float sum = 0;
1955
1956             for (int k = 0; k < 3; k++)
1957                 sum += a[i][k] * b[k][j];
1958
1959             c[i][j] = sum;
1960         }
1961     }
1962 }
1963
1964 /**
1965  * Calculate rotation matrix for yaw/pitch/roll angles.
1966  */
1967 static inline void calculate_rotation_matrix(float yaw, float pitch, float roll,
1968                                              float rot_mat[3][3],
1969                                              const int rotation_order[3])
1970 {
1971     const float yaw_rad   = yaw   * M_PI / 180.f;
1972     const float pitch_rad = pitch * M_PI / 180.f;
1973     const float roll_rad  = roll  * M_PI / 180.f;
1974
1975     const float sin_yaw   = sinf(-yaw_rad);
1976     const float cos_yaw   = cosf(-yaw_rad);
1977     const float sin_pitch = sinf(pitch_rad);
1978     const float cos_pitch = cosf(pitch_rad);
1979     const float sin_roll  = sinf(roll_rad);
1980     const float cos_roll  = cosf(roll_rad);
1981
1982     float m[3][3][3];
1983     float temp[3][3];
1984
1985     m[0][0][0] =  cos_yaw;  m[0][0][1] = 0;          m[0][0][2] =  sin_yaw;
1986     m[0][1][0] =  0;        m[0][1][1] = 1;          m[0][1][2] =  0;
1987     m[0][2][0] = -sin_yaw;  m[0][2][1] = 0;          m[0][2][2] =  cos_yaw;
1988
1989     m[1][0][0] = 1;         m[1][0][1] = 0;          m[1][0][2] =  0;
1990     m[1][1][0] = 0;         m[1][1][1] = cos_pitch;  m[1][1][2] = -sin_pitch;
1991     m[1][2][0] = 0;         m[1][2][1] = sin_pitch;  m[1][2][2] =  cos_pitch;
1992
1993     m[2][0][0] = cos_roll;  m[2][0][1] = -sin_roll;  m[2][0][2] =  0;
1994     m[2][1][0] = sin_roll;  m[2][1][1] =  cos_roll;  m[2][1][2] =  0;
1995     m[2][2][0] = 0;         m[2][2][1] =  0;         m[2][2][2] =  1;
1996
1997     multiply_matrix(temp, m[rotation_order[0]], m[rotation_order[1]]);
1998     multiply_matrix(rot_mat, temp, m[rotation_order[2]]);
1999 }
2000
2001 /**
2002  * Rotate vector with given rotation matrix.
2003  *
2004  * @param rot_mat rotation matrix
2005  * @param vec vector
2006  */
2007 static inline void rotate(const float rot_mat[3][3],
2008                           float *vec)
2009 {
2010     const float x_tmp = vec[0] * rot_mat[0][0] + vec[1] * rot_mat[0][1] + vec[2] * rot_mat[0][2];
2011     const float y_tmp = vec[0] * rot_mat[1][0] + vec[1] * rot_mat[1][1] + vec[2] * rot_mat[1][2];
2012     const float z_tmp = vec[0] * rot_mat[2][0] + vec[1] * rot_mat[2][1] + vec[2] * rot_mat[2][2];
2013
2014     vec[0] = x_tmp;
2015     vec[1] = y_tmp;
2016     vec[2] = z_tmp;
2017 }
2018
2019 static inline void set_mirror_modifier(int h_flip, int v_flip, int d_flip,
2020                                        float *modifier)
2021 {
2022     modifier[0] = h_flip ? -1.f : 1.f;
2023     modifier[1] = v_flip ? -1.f : 1.f;
2024     modifier[2] = d_flip ? -1.f : 1.f;
2025 }
2026
2027 static inline void mirror(const float *modifier, float *vec)
2028 {
2029     vec[0] *= modifier[0];
2030     vec[1] *= modifier[1];
2031     vec[2] *= modifier[2];
2032 }
2033
2034 static int allocate_plane(V360Context *s, int sizeof_uv, int sizeof_ker, int p)
2035 {
2036     s->u[p] = av_calloc(s->planewidth[p] * s->planeheight[p], sizeof_uv);
2037     s->v[p] = av_calloc(s->planewidth[p] * s->planeheight[p], sizeof_uv);
2038     if (!s->u[p] || !s->v[p])
2039         return AVERROR(ENOMEM);
2040     if (sizeof_ker) {
2041         s->ker[p] = av_calloc(s->planewidth[p] * s->planeheight[p], sizeof_ker);
2042         if (!s->ker[p])
2043             return AVERROR(ENOMEM);
2044     }
2045
2046     return 0;
2047 }
2048
2049 static int config_output(AVFilterLink *outlink)
2050 {
2051     AVFilterContext *ctx = outlink->src;
2052     AVFilterLink *inlink = ctx->inputs[0];
2053     V360Context *s = ctx->priv;
2054     const AVPixFmtDescriptor *desc = av_pix_fmt_desc_get(inlink->format);
2055     const int depth = desc->comp[0].depth;
2056     int sizeof_uv;
2057     int sizeof_ker;
2058     int elements;
2059     int err;
2060     int p, h, w;
2061     float hf, wf;
2062     float output_mirror_modifier[3];
2063     void (*in_transform)(const V360Context *s,
2064                          const float *vec, int width, int height,
2065                          uint16_t us[4][4], uint16_t vs[4][4], float *du, float *dv);
2066     void (*out_transform)(const V360Context *s,
2067                           int i, int j, int width, int height,
2068                           float *vec);
2069     void (*calculate_kernel)(float du, float dv, const XYRemap *r_tmp,
2070                              uint16_t *u, uint16_t *v, int16_t *ker);
2071     float rot_mat[3][3];
2072
2073     s->input_mirror_modifier[0] = s->ih_flip ? -1.f : 1.f;
2074     s->input_mirror_modifier[1] = s->iv_flip ? -1.f : 1.f;
2075
2076     switch (s->interp) {
2077     case NEAREST:
2078         calculate_kernel = nearest_kernel;
2079         s->remap_slice = depth <= 8 ? remap1_8bit_slice : remap1_16bit_slice;
2080         elements = 1;
2081         sizeof_uv = sizeof(uint16_t) * elements;
2082         sizeof_ker = 0;
2083         break;
2084     case BILINEAR:
2085         calculate_kernel = bilinear_kernel;
2086         s->remap_slice = depth <= 8 ? remap2_8bit_slice : remap2_16bit_slice;
2087         elements = 2 * 2;
2088         sizeof_uv = sizeof(uint16_t) * elements;
2089         sizeof_ker = sizeof(uint16_t) * elements;
2090         break;
2091     case BICUBIC:
2092         calculate_kernel = bicubic_kernel;
2093         s->remap_slice = depth <= 8 ? remap4_8bit_slice : remap4_16bit_slice;
2094         elements = 4 * 4;
2095         sizeof_uv = sizeof(uint16_t) * elements;
2096         sizeof_ker = sizeof(uint16_t) * elements;
2097         break;
2098     case LANCZOS:
2099         calculate_kernel = lanczos_kernel;
2100         s->remap_slice = depth <= 8 ? remap4_8bit_slice : remap4_16bit_slice;
2101         elements = 4 * 4;
2102         sizeof_uv = sizeof(uint16_t) * elements;
2103         sizeof_ker = sizeof(uint16_t) * elements;
2104         break;
2105     default:
2106         av_assert0(0);
2107     }
2108
2109     ff_v360_init(s, depth);
2110
2111     for (int order = 0; order < NB_RORDERS; order++) {
2112         const char c = s->rorder[order];
2113         int rorder;
2114
2115         if (c == '\0') {
2116             av_log(ctx, AV_LOG_ERROR,
2117                    "Incomplete rorder option. Direction for all 3 rotation orders should be specified.\n");
2118             return AVERROR(EINVAL);
2119         }
2120
2121         rorder = get_rorder(c);
2122         if (rorder == -1) {
2123             av_log(ctx, AV_LOG_ERROR,
2124                    "Incorrect rotation order symbol '%c' in rorder option.\n", c);
2125             return AVERROR(EINVAL);
2126         }
2127
2128         s->rotation_order[order] = rorder;
2129     }
2130
2131     switch (s->in) {
2132     case EQUIRECTANGULAR:
2133         in_transform = xyz_to_equirect;
2134         err = 0;
2135         wf = inlink->w;
2136         hf = inlink->h;
2137         break;
2138     case CUBEMAP_3_2:
2139         in_transform = xyz_to_cube3x2;
2140         err = prepare_cube_in(ctx);
2141         wf = inlink->w / 3.f * 4.f;
2142         hf = inlink->h;
2143         break;
2144     case CUBEMAP_1_6:
2145         in_transform = xyz_to_cube1x6;
2146         err = prepare_cube_in(ctx);
2147         wf = inlink->w * 4.f;
2148         hf = inlink->h / 3.f;
2149         break;
2150     case CUBEMAP_6_1:
2151         in_transform = xyz_to_cube6x1;
2152         err = prepare_cube_in(ctx);
2153         wf = inlink->w / 3.f * 2.f;
2154         hf = inlink->h * 2.f;
2155         break;
2156     case EQUIANGULAR:
2157         in_transform = xyz_to_eac;
2158         err = prepare_eac_in(ctx);
2159         wf = inlink->w;
2160         hf = inlink->h / 9.f * 8.f;
2161         break;
2162     case FLAT:
2163         av_log(ctx, AV_LOG_ERROR, "Flat format is not accepted as input.\n");
2164         return AVERROR(EINVAL);
2165     case DUAL_FISHEYE:
2166         in_transform = xyz_to_dfisheye;
2167         err = 0;
2168         wf = inlink->w;
2169         hf = inlink->h;
2170         break;
2171     case BARREL:
2172         in_transform = xyz_to_barrel;
2173         err = 0;
2174         wf = inlink->w / 5.f * 4.f;
2175         hf = inlink->h;
2176         break;
2177     default:
2178         av_log(ctx, AV_LOG_ERROR, "Specified input format is not handled.\n");
2179         return AVERROR_BUG;
2180     }
2181
2182     if (err != 0) {
2183         return err;
2184     }
2185
2186     switch (s->out) {
2187     case EQUIRECTANGULAR:
2188         out_transform = equirect_to_xyz;
2189         err = 0;
2190         w = roundf(wf);
2191         h = roundf(hf);
2192         break;
2193     case CUBEMAP_3_2:
2194         out_transform = cube3x2_to_xyz;
2195         err = prepare_cube_out(ctx);
2196         w = roundf(wf / 4.f * 3.f);
2197         h = roundf(hf);
2198         break;
2199     case CUBEMAP_1_6:
2200         out_transform = cube1x6_to_xyz;
2201         err = prepare_cube_out(ctx);
2202         w = roundf(wf / 4.f);
2203         h = roundf(hf * 3.f);
2204         break;
2205     case CUBEMAP_6_1:
2206         out_transform = cube6x1_to_xyz;
2207         err = prepare_cube_out(ctx);
2208         w = roundf(wf / 2.f * 3.f);
2209         h = roundf(hf / 2.f);
2210         break;
2211     case EQUIANGULAR:
2212         out_transform = eac_to_xyz;
2213         err = prepare_eac_out(ctx);
2214         w = roundf(wf);
2215         h = roundf(hf / 8.f * 9.f);
2216         break;
2217     case FLAT:
2218         out_transform = flat_to_xyz;
2219         err = prepare_flat_out(ctx);
2220         w = roundf(wf * s->flat_range[0] / s->flat_range[1] / 2.f);
2221         h = roundf(hf);
2222         break;
2223     case DUAL_FISHEYE:
2224         av_log(ctx, AV_LOG_ERROR, "Dual fisheye format is not accepted as output.\n");
2225         return AVERROR(EINVAL);
2226     case BARREL:
2227         out_transform = barrel_to_xyz;
2228         err = 0;
2229         w = roundf(wf / 4.f * 5.f);
2230         h = roundf(hf);
2231         break;
2232     case STEREOGRAPHIC:
2233         out_transform = stereographic_to_xyz;
2234         err = 0;
2235         w = FFMAX(roundf(wf), roundf(hf));
2236         h = w;
2237         break;
2238     default:
2239         av_log(ctx, AV_LOG_ERROR, "Specified output format is not handled.\n");
2240         return AVERROR_BUG;
2241     }
2242
2243     if (err != 0) {
2244         return err;
2245     }
2246
2247     // Override resolution with user values if specified
2248     if (s->width > 0 && s->height > 0) {
2249         w = s->width;
2250         h = s->height;
2251     } else if (s->width > 0 || s->height > 0) {
2252         av_log(ctx, AV_LOG_ERROR, "Both width and height values should be specified.\n");
2253         return AVERROR(EINVAL);
2254     } else {
2255         if (s->out_transpose)
2256             FFSWAP(int, w, h);
2257
2258         if (s->in_transpose)
2259             FFSWAP(int, w, h);
2260     }
2261
2262     s->planeheight[1] = s->planeheight[2] = FF_CEIL_RSHIFT(h, desc->log2_chroma_h);
2263     s->planeheight[0] = s->planeheight[3] = h;
2264     s->planewidth[1]  = s->planewidth[2] = FF_CEIL_RSHIFT(w, desc->log2_chroma_w);
2265     s->planewidth[0]  = s->planewidth[3] = w;
2266
2267     outlink->h = h;
2268     outlink->w = w;
2269
2270     s->inplaneheight[1] = s->inplaneheight[2] = FF_CEIL_RSHIFT(inlink->h, desc->log2_chroma_h);
2271     s->inplaneheight[0] = s->inplaneheight[3] = inlink->h;
2272     s->inplanewidth[1]  = s->inplanewidth[2]  = FF_CEIL_RSHIFT(inlink->w, desc->log2_chroma_w);
2273     s->inplanewidth[0]  = s->inplanewidth[3]  = inlink->w;
2274     s->nb_planes = av_pix_fmt_count_planes(inlink->format);
2275
2276     if (desc->log2_chroma_h == desc->log2_chroma_w && desc->log2_chroma_h == 0) {
2277         s->nb_allocated = 1;
2278         s->map[0] = s->map[1] = s->map[2] = s->map[3] = 0;
2279         allocate_plane(s, sizeof_uv, sizeof_ker, 0);
2280     } else if (desc->log2_chroma_h == desc->log2_chroma_w) {
2281         s->nb_allocated = 2;
2282         s->map[0] = 0;
2283         s->map[1] = s->map[2] = 1;
2284         s->map[3] = 0;
2285         allocate_plane(s, sizeof_uv, sizeof_ker, 0);
2286         allocate_plane(s, sizeof_uv, sizeof_ker, 1);
2287     } else {
2288         s->nb_allocated = 3;
2289         s->map[0] = 0;
2290         s->map[1] = 1;
2291         s->map[2] = 2;
2292         s->map[3] = 0;
2293         allocate_plane(s, sizeof_uv, sizeof_ker, 0);
2294         allocate_plane(s, sizeof_uv, sizeof_ker, 1);
2295         allocate_plane(s, sizeof_uv, sizeof_ker, 2);
2296     }
2297
2298     calculate_rotation_matrix(s->yaw, s->pitch, s->roll, rot_mat, s->rotation_order);
2299     set_mirror_modifier(s->h_flip, s->v_flip, s->d_flip, output_mirror_modifier);
2300
2301     // Calculate remap data
2302     for (p = 0; p < s->nb_allocated; p++) {
2303         const int width = s->planewidth[p];
2304         const int height = s->planeheight[p];
2305         const int in_width = s->inplanewidth[p];
2306         const int in_height = s->inplaneheight[p];
2307         float du, dv;
2308         float vec[3];
2309         XYRemap r_tmp;
2310         int i, j;
2311
2312         for (i = 0; i < width; i++) {
2313             for (j = 0; j < height; j++) {
2314                 uint16_t *u = s->u[p] + (j * width + i) * elements;
2315                 uint16_t *v = s->v[p] + (j * width + i) * elements;
2316                 int16_t *ker = s->ker[p] + (j * width + i) * elements;
2317
2318                 if (s->out_transpose)
2319                     out_transform(s, j, i, height, width, vec);
2320                 else
2321                     out_transform(s, i, j, width, height, vec);
2322                 rotate(rot_mat, vec);
2323                 mirror(output_mirror_modifier, vec);
2324                 if (s->in_transpose)
2325                     in_transform(s, vec, in_height, in_width, r_tmp.v, r_tmp.u, &du, &dv);
2326                 else
2327                     in_transform(s, vec, in_width, in_height, r_tmp.u, r_tmp.v, &du, &dv);
2328                 calculate_kernel(du, dv, &r_tmp, u, v, ker);
2329             }
2330         }
2331     }
2332
2333     return 0;
2334 }
2335
2336 static int filter_frame(AVFilterLink *inlink, AVFrame *in)
2337 {
2338     AVFilterContext *ctx = inlink->dst;
2339     AVFilterLink *outlink = ctx->outputs[0];
2340     V360Context *s = ctx->priv;
2341     AVFrame *out;
2342     ThreadData td;
2343
2344     out = ff_get_video_buffer(outlink, outlink->w, outlink->h);
2345     if (!out) {
2346         av_frame_free(&in);
2347         return AVERROR(ENOMEM);
2348     }
2349     av_frame_copy_props(out, in);
2350
2351     td.in = in;
2352     td.out = out;
2353
2354     ctx->internal->execute(ctx, s->remap_slice, &td, NULL, FFMIN(outlink->h, ff_filter_get_nb_threads(ctx)));
2355
2356     av_frame_free(&in);
2357     return ff_filter_frame(outlink, out);
2358 }
2359
2360 static av_cold void uninit(AVFilterContext *ctx)
2361 {
2362     V360Context *s = ctx->priv;
2363     int p;
2364
2365     for (p = 0; p < s->nb_allocated; p++) {
2366         av_freep(&s->u[p]);
2367         av_freep(&s->v[p]);
2368         av_freep(&s->ker[p]);
2369     }
2370 }
2371
2372 static const AVFilterPad inputs[] = {
2373     {
2374         .name         = "default",
2375         .type         = AVMEDIA_TYPE_VIDEO,
2376         .filter_frame = filter_frame,
2377     },
2378     { NULL }
2379 };
2380
2381 static const AVFilterPad outputs[] = {
2382     {
2383         .name         = "default",
2384         .type         = AVMEDIA_TYPE_VIDEO,
2385         .config_props = config_output,
2386     },
2387     { NULL }
2388 };
2389
2390 AVFilter ff_vf_v360 = {
2391     .name          = "v360",
2392     .description   = NULL_IF_CONFIG_SMALL("Convert 360 projection of video."),
2393     .priv_size     = sizeof(V360Context),
2394     .uninit        = uninit,
2395     .query_formats = query_formats,
2396     .inputs        = inputs,
2397     .outputs       = outputs,
2398     .priv_class    = &v360_class,
2399     .flags         = AVFILTER_FLAG_SLICE_THREADS,
2400 };