#@gmic scale_dcci2x : _edge_threshold>=0
#@gmic : Double image size using directional cubic convolution interpolation,
#@gmic : as described in https://en.wikipedia.org/wiki/Directional_Cubic_Convolution_Interpolation.
#@gmic : Default value: 'edge_threshold=1.15'.
#@gmic : $ image.jpg --scale_dcci2x ,
scale_dcci2x : -check "${1=1.15}>=0"
  -e[^-1] "Double image size using directional cubic convolution interpolation."
  -v - -repeat $! -l[$>]
    -r 200%,200%,1,100%,4 -r {w-1},{h-1},1,100%,0
    -f "*
if (!(x%2 && y%2),i,
  d1 = abs(j(-1,-3,0,0,0,1) - j(-3,-1,0,0,0,1)) + abs(j(1,-3,0,0,0,1) - j(-1,-1,0,0,0,1)) + abs(j(3,-3,0,0,0,1) - j(1,-1,0,0,0,1)) +
       abs(j(-1,-1,0,0,0,1) - j(-3,1,0,0,0,1)) + abs(j(1,-1,0,0,0,1) - j(-1,1,0,0,0,1)) + abs(j(3,-1,0,0,0,1) - j(1,1,0,0,0,1)) +
       abs(j(-1,1,0,0,0,1) - j(-3,3,0,0,0,1)) + abs(j(1,1,0,0,0,1) - j(-1,3,0,0,0,1)) + abs(j(3,1,0,0,0,1) - j(1,3,0,0,0,1));
  d2 = abs(j(-3,-3,0,0,0,1) - j(-1,-1,0,0,0,1)) + abs(j(-1,-3,0,0,0,1) - j(1,-1,0,0,0,1)) + abs(j(1,-3,0,0,0,1) - j(3,-1,0,0,0,1)) +
       abs(j(-3,-1,0,0,0,1) - j(-1,1,0,0,0,1)) + abs(j(-1,-1,0,0,0,1) - j(1,1,0,0,0,1)) + abs(j(1,-1,0,0,0,1) - j(3,1,0,0,0,1)) +
       abs(j(-3,1,0,0,0,1) - j(-1,3,0,0,0,1)) + abs(j(-1,1,0,0,0,1) - j(1,3,0,0,0,1)) + abs(j(1,1,0,0,0,1) - j(3,3,0,0,0,1));
  if ((1 + d1)/(1 + d2)>$1,
    (-j(-3,-3,0,0,0,1) + 9*j(-1,-1,0,0,0,1) + 9*j(1,1,0,0,0,1) - j(3,3,0,0,0,1))/16,
  if ((1 + d2)/(1 + d1)>$1,
    (-j(3,-3,0,0,0,1) + 9*j(1,-1,0,0,0,1) + 9*j(-1,1,0,0,0,1) - j(-3,3,0,0,0,1))/16,
    (
      w1 = 1/(1 + d1^5);
      w2 = 1/(1 + d2^5);
      weight1 = w1/(w1 + w2);
      weight2 = w2/(w1 + w2);
      DownRightPixel = (-j(-3,-3,0,0,0,1) + 9*j(-1,-1,0,0,0,1) + 9*j(1,1,0,0,0,1) - j(3,3,0,0,0,1))/16;
      UpRightPixel = (-j(3,-3,0,0,0,1) + 9*j(1,-1,0,0,0,1) + 9*j(-1,1,0,0,0,1) - j(-3,3,0,0,0,1))/16;
      DownRightPixel*weight1 + UpRightPixel*weight2
    )
  ))
)
"
    -f "*
if ((!(x%2) && !(y%2)) || ((x%2) && (y%2)),i,
  d1 = abs(j(1,-2,0,0,0,1) - j(-1,-2,0,0,0,1)) +
       abs(j(2,-1,0,0,0,1) - j(0,-1,0,0,0,1)) + abs(j(0,-1,0,0,0,1) - j(-2,-1,0,0,0,1)) +
       abs(j(3,0,0,0,0,1) - j(1,0,0,0,0,1)) + abs(j(1,0,0,0,0,1) - j(-1,0,0,0,0,1)) + abs(j(-1,0,0,0,0,1) - j(-3,0,0,0,0,1)) +
       abs(j(2,1,0,0,0,1) - j(0,1,0,0,0,1)) + abs(j(0,1,0,0,0,1) - j(-2,1,0,0,0,1)) +
       abs(j(1,2,0,0,0,1) - j(-1,2,0,0,0,1));
  d2 = abs(j(-2,1,0,0,0,1) - j(-2,-1,0,0,0,1)) +
       abs(j(-1,2,0,0,0,1) - j(-1,0,0,0,0,1)) + abs(j(-1,0,0,0,0,1) - j(-1,-2,0,0,0,1)) +
       abs(j(0,3,0,0,0,1) - j(0,1,0,0,0,1)) + abs(j(0,1,0,0,0,1) - j(0,-1,0,0,0,1)) + abs(j(0,-1,0,0,0,1) - j(0,-3,0,0,0,1)) +
       abs(j(1,2,0,0,0,1) - j(1,0,0,0,0,1)) + abs(j(1,0,0,0,0,1) - j(1,-2,0,0,0,1)) +
       abs(j(2,1,0,0,0,1) - j(2,-1,0,0,0,1));
  if ((1 + d1)/(1 + d2)>$1,
    (-j(0,-3,0,0,0,1) + 9*j(0,-1,0,0,0,1) + 9*j(0,1,0,0,0,1) - j(0,3,0,0,0,1))/16,
  if ((1 + d2)/(1 + d1)>$1,
    (-j(-3,0,0,0,0,1) + 9*j(-1,0,0,0,0,1) + 9*j(1,0,0,0,0,1) - j(3,0,0,0,0,1))/16,
    (
      w1 = 1/(1 + d1^5);
      w2 = 1/(1 + d2^5);
      weight1 = w1/(w1 + w2);
      weight2 = w2/(w1 + w2);
      HorizontalPixel = (-j(-3,0,0,0,0,1) + 9*j(-1,0,0,0,0,1) + 9*j(1,0,0,0,0,1) - j(3,0,0,0,0,1))/16;
      VerticalPixel = (-j(0,-3,0,0,0,1) + 9*j(0,-1,0,0,0,1) + 9*j(0,1,0,0,0,1) - j(0,3,0,0,0,1))/16;
      VerticalPixel*weight1 + HorizontalPixel*weight2
    )
  ))
)
"
  -endl -done -v +