{
 "cells": [
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "## Vertical stripe fix for 5D Mark III H.264 files"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "### What is this about?\n",
    "\n",
    "Some users noticed a vertical pattern in 5D Mark III video files, first with regular raw video images. Some converters will correct this artifact automatically, so it wasn't a big issue lately.\n",
    "\n",
    "After implementing the [3x crop mode](http://magiclantern.fm/forum/index.php?topic=17021), the artifact became visible in H.264 files as well. The previous algorithm that operates on Bayer raw data can no longer be used - it need to be adapted to H.264 files.\n",
    "\n",
    "The effect is visible only in highlights, not in shadows. If you see vertical lines visible in shadows, that's a different issue, so you may stop reading here.\n",
    "\n",
    "### What's causing it?\n",
    "\n",
    "The 5D Mark III appears to read out 8 columns at a time, in parallel. These 8 columns appear to have different amplifier gains, and that's why the effect is visible in highlights, but not in shadows.\n",
    "\n",
    "Since:\n",
    "\n",
    "1. with regular (non-crop) mode, the artifact was visible in RAW, but not in H.264, \n",
    "2. with crop mode, the artifact is visible in H.264, but not in RAW,\n",
    "3. the crop mode only changes the sampled area on the sensor, and the pixel binning factors, but leaves all other parameters (resolution, Canon processing and so on) unchanged,\n",
    "\n",
    "we may suspect that those stripes are present in the raw sensor data, corrected by Canon code in regular (non-crop) H.264, but for crop mode H.264, it would require a recalibration (and reverse engineering to figure out how to do that first). It is also possible that our vertical stripe artifact is amplified by the mismatched calibration.\n",
    "\n",
    "\n",
    "### Enough chit-chat, how to fix it?\n",
    "\n",
    "You will need:\n",
    "\n",
    "- IPython notebook\n",
    "- octave (the one prepackaged for your operating system should be fine)\n",
    "- the 'image' package from octave-forge\n",
    "- octave_kernel for the IPython notebook\n",
    "- ffmpeg\n",
    "\n",
    "Let's the sample files from [here](http://magiclantern.fm/forum/index.php?topic=17021.msg166405#msg166405) in octave. We'll extract a single frame from each clip, as uncompressed 8-bit ppm:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%%shell\n",
    "# note: ffmpeg is very verbose, so I won't display the output here\n",
    "ffmpeg -ss 0.96 -i 1080p\\ Stripes\\ 1.mov bad.ppm -y 2>ffmpeg.log\n",
    "ffmpeg -i 1080p\\ Stripes\\ 2.mov ref.ppm -y 2>>ffmpeg.log"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Octave code follows:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "pkg load image\n",
    "more off"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "bad = im2double(imread('bad.ppm'));\n",
    "ref = im2double(imread('ref.ppm'));\n",
    "figure,imshow(bad(1:5:end,1:5:end,:))\n",
    "figure,imshow(ref(1:5:end,1:5:end,:))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The user who reported the issue included a normal image and also a blank wall, that can be used to extract the vertical pattern. Nice!\n",
    "\n",
    "Let's see a 1:1 crop from the bad image:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "imshow(bad(1:400,1:800,:))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "Doesn't look very nice.\n",
    "\n",
    "Let's try to learn the pattern from the reference image (blank wall), from the green channel. Easier to work with grayscale, isn't?"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "% will stretch the image to make the pattern very obvious\n",
    "ref_g = ref(:,:,2);\n",
    "imshow(ref_g(1:400,1:800), [])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "The reference blank wall is not completely uniform - let's filter out the low frequency components:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "% filter out very low frequency components in the horizontal direction\n",
    "im = ref_g;\n",
    "imf = imfilter(im, ones(1,100)/100, 'replicate');\n",
    "figure, imshow((im-imf)(1:400,1:800), [])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {
    "collapsed": true
   },
   "source": [
    "The image is more uniform now, but there are still a bunch of artifacts. Let's find out the gain for each column:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "% get the vertical pattern\n",
    "% we already know the pattern is actually a difference in column gains\n",
    "% assume the filtered image is ideal, and extract per-column gain variations\n",
    "gain = median(im ./ imf);\n",
    "plot(gain)\n",
    "axis tight"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's put it together:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "function cor = identify_pattern_gray(im)\n",
    "  % filter out very low frequency components in the horizontal direction\n",
    "  imf = imfilter(im, ones(1,100)/100, 'replicate');\n",
    "\n",
    "  % get the vertical pattern\n",
    "  % we already know the pattern is actually a difference in column gains\n",
    "  % assume the filtered image is ideal, and extract per-column gain variations\n",
    "  gain = median(im ./ imf);\n",
    "  \n",
    "  % expand the correction data (vector, one entry per column)\n",
    "  % to match the size of the image\n",
    "  cor = bsxfun(@plus, im*0, gain);\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "cor_g = identify_pattern_gray(ref_g);\n",
    "imshow(cor_g(1:400,1:800), []);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That is the correction pattern - stretched to see it better.\n",
    "\n",
    "You probably noticed it doesn't repeat every 8 columns, as you would expect from my initial description. In raw, it actually does. However, the H.264 1920 appears to be upsampled from 1904, and this explains the not-exactly-periodic pattern you have just seen.\n",
    "\n",
    "The above image (and all the other crops from this page) are 1:1, non-resized. What you are seeing **is** the actual pattern that affects the image.\n",
    "\n",
    "Let's check the corrected green channel:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "figure, imshow((ref_g ./ cor_g)(1:400,1:800), []);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks like it worked!\n",
    "\n",
    "What are those ugly blocks? H.264 compression artifacts :)\n",
    "\n",
    "Let's mask them with a bit of noise:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    " imshow((ref_g ./ cor_g + randn(size(ref_g))/200)(1:400,1:800), [])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Which one do you prefer?\n",
    "\n",
    "-----------------"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Forget the noise for now; let's check the other color channels:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "ref_r = ref(:,:,1);\n",
    "ref_b = ref(:,:,3);\n",
    "cor_r = identify_pattern_gray(ref_r, 0);\n",
    "cor_b = identify_pattern_gray(ref_b, 0);\n",
    "figure,imshow([cor_r(1:100,1:800); cor_g(1:100,1:800); cor_b(1:100,1:800)], []);"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "So it looks like every channel needs its own correction. No big deal.\n",
    "\n",
    "Let's put it together."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "function cor = identify_pattern(im)\n",
    "  for c = 1:size(im,3)\n",
    "    cor(:,:,c) = identify_pattern_gray(im(:,:,c));\n",
    "  end\n",
    "end"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "cor = identify_pattern(ref);\n",
    "fix = bad ./ cor;\n",
    "\n",
    "imshow(fix(1:400,1:800,:))"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Not bad :)\n",
    "\n",
    "Let's pixel-peep the blank wall a little more, on the green channel:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "imshow(fix(1:200,501:1000,2), [])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Let's add the noise again:"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    " fixn = fix + randn(size(fix))/200;\n",
    " figure,imshow(fixn(1:400,1:800,:))\n",
    " figure,imshow(fixn(1:200,501:1000,2), [min(fix(1:200,501:1000,2)(:)) max(fix(1:200,501:1000,2)(:))])"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Looks fine. Now let's save the two results to files."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": true
   },
   "outputs": [],
   "source": [
    "imwrite(fix, 'fix.jpg', 'quality', 99)\n",
    "imwrite(fixn, 'fixn.jpg', 'quality', 99)"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Now let's process all the frames from our short sample clip. I like the noisy version better, so I'll just use that."
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%%shell\n",
    "ffmpeg -i 1080p\\ Stripes\\ 1.mov frame%03d.ppm -y 2>ffmpeg.log"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "for i = 1:25\n",
    "   fin = sprintf('frame%03d.ppm', i)\n",
    "   fout = sprintf('fixed%03d.ppm', i);\n",
    "   im = im2double(imread(fin));\n",
    "   fix = im ./ cor + randn(size(im))/200;\n",
    "   imwrite(fix, fout);\n",
    "end"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "Assemble the frames (feel free to change the codec, I'm not familiar with them):"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false
   },
   "outputs": [],
   "source": [
    "%%shell\n",
    "ffmpeg -i fixed%03d.ppm -pix_fmt yuv420p -vcodec rawvideo video.avi -y 2>ffmpeg.log"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "And now let's admire our masterpiece :)"
   ]
  },
  {
   "cell_type": "code",
   "execution_count": null,
   "metadata": {
    "collapsed": false,
    "scrolled": true
   },
   "outputs": [],
   "source": [
    "%%shell\n",
    "ffplay -loop 100 video-n.avi 2>ffmpeg.log"
   ]
  },
  {
   "cell_type": "markdown",
   "metadata": {},
   "source": [
    "That's all, folks!\n",
    "\n",
    "Hope this notebook contains all you need to integrate the method in your video workflow, and to try new ideas."
   ]
  }
 ],
 "metadata": {
  "kernelspec": {
   "display_name": "Octave",
   "language": "octave",
   "name": "octave"
  },
  "language_info": {
   "file_extension": ".m",
   "help_links": [
    {
     "text": "MetaKernel Magics",
     "url": "https://github.com/calysto/metakernel/blob/master/metakernel/magics/README.md"
    }
   ],
   "mimetype": "text/x-octave",
   "name": "octave",
   "version": "0.15.10"
  }
 },
 "nbformat": 4,
 "nbformat_minor": 0
}