{
  "cells": [
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "\n# Image resampling\n\nImages are represented by discrete pixels assigned color values, either on the\nscreen or in an image file.  When a user calls `~.Axes.imshow` with a data\narray, it is rare that the size of the data array exactly matches the number of\npixels allotted to the image in the figure, so Matplotlib resamples or [scales](https://en.wikipedia.org/wiki/Image_scaling) the data or image to fit.  If\nthe data array is larger than the number of pixels allotted in the rendered figure,\nthen the image will be \"down-sampled\" and image information will be lost.\nConversely, if the data array is smaller than the number of output pixels then each\ndata point will get multiple pixels, and the image is \"up-sampled\".\n\nIn the following figure, the first data array has size (450, 450), but is\nrepresented by far fewer pixels in the figure, and hence is down-sampled.  The\nsecond data array has size (4, 4), and is represented by far more pixels, and\nhence is up-sampled.\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "import matplotlib.pyplot as plt\nimport numpy as np\n\nfig, axs = plt.subplots(1, 2, figsize=(4, 2))\n\n# First we generate a 450x450 pixel image with varying frequency content:\nN = 450\nx = np.arange(N) / N - 0.5\ny = np.arange(N) / N - 0.5\naa = np.ones((N, N))\naa[::2, :] = -1\n\nX, Y = np.meshgrid(x, y)\nR = np.sqrt(X**2 + Y**2)\nf0 = 5\nk = 100\na = np.sin(np.pi * 2 * (f0 * R + k * R**2 / 2))\n# make the left hand side of this\na[:int(N / 2), :][R[:int(N / 2), :] < 0.4] = -1\na[:int(N / 2), :][R[:int(N / 2), :] < 0.3] = 1\naa[:, int(N / 3):] = a[:, int(N / 3):]\nalarge = aa\n\naxs[0].imshow(alarge, cmap='RdBu_r')\naxs[0].set_title('(450, 450) Down-sampled', fontsize='medium')\n\nnp.random.seed(19680801+9)\nasmall = np.random.rand(4, 4)\naxs[1].imshow(asmall, cmap='viridis')\naxs[1].set_title('(4, 4) Up-sampled', fontsize='medium')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Matplotlib's `~.Axes.imshow` method has two keyword arguments to allow the user\nto control how resampling is done.  The *interpolation* keyword argument allows\na choice of the kernel that is used for resampling, allowing either [anti-alias](https://en.wikipedia.org/wiki/Anti-aliasing_filter) filtering if\ndown-sampling, or smoothing of pixels if up-sampling.  The\n*interpolation_stage* keyword argument, determines if this smoothing kernel is\napplied to the underlying data, or if the kernel is applied to the RGBA pixels.\n\n``interpolation_stage='rgba'``: Data -> Normalize -> RGBA -> Interpolate/Resample\n\n``interpolation_stage='data'``: Data -> Interpolate/Resample -> Normalize -> RGBA\n\nFor both keyword arguments, Matplotlib has a default \"auto\", that is\nrecommended for most situations, and is described below.  Note that this\ndefault behaves differently if the image is being down- or up-sampled, as\ndescribed below.\n\n## Down-sampling and modest up-sampling\n\nWhen down-sampling data, we usually want to remove aliasing by smoothing the\nimage first and then sub-sampling it.  In Matplotlib, we can do that smoothing\nbefore mapping the data to colors, or we can do the smoothing on the RGB(A)\nimage pixels.  The differences between these are shown below, and controlled\nwith the *interpolation_stage* keyword argument.\n\nThe following images are down-sampled from 450 data pixels to approximately\n125 pixels or 250 pixels (depending on your display).\nThe underlying image has alternating +1, -1 stripes on the left side, and\na varying wavelength ([chirp](https://en.wikipedia.org/wiki/Chirp)) pattern\nin the rest of the image.  If we zoom, we can see this detail without any\ndown-sampling:\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fig, ax = plt.subplots(figsize=(4, 4), layout='compressed')\nax.imshow(alarge, interpolation='nearest', cmap='RdBu_r')\nax.set_xlim(100, 200)\nax.set_ylim(275, 175)\nax.set_title('Zoom')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "If we down-sample, the simplest algorithm is to decimate the data using\n[nearest-neighbor interpolation](https://en.wikipedia.org/wiki/Nearest-neighbor_interpolation).  We can\ndo this in either data space or RGBA space:\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fig, axs = plt.subplots(1, 2, figsize=(5, 2.7), layout='compressed')\nfor ax, interp, space in zip(axs.flat, ['nearest', 'nearest'],\n                                       ['data', 'rgba']):\n    ax.imshow(alarge, interpolation=interp, interpolation_stage=space,\n              cmap='RdBu_r')\n    ax.set_title(f\"interpolation='{interp}'\\nstage='{space}'\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Nearest interpolation is identical in data and RGBA space, and both exhibit\n[Moir\u00e9](https://en.wikipedia.org/wiki/Moir\u00e9_pattern) patterns because the\nhigh-frequency data is being down-sampled and shows up as lower frequency\npatterns. We can reduce the Moir\u00e9 patterns by applying an anti-aliasing filter\nto the image before rendering:\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fig, axs = plt.subplots(1, 2, figsize=(5, 2.7), layout='compressed')\nfor ax, interp, space in zip(axs.flat, ['hanning', 'hanning'],\n                                       ['data', 'rgba']):\n    ax.imshow(alarge, interpolation=interp, interpolation_stage=space,\n              cmap='RdBu_r')\n    ax.set_title(f\"interpolation='{interp}'\\nstage='{space}'\")\nplt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The [Hanning](https://en.wikipedia.org/wiki/Hann_function) filter smooths\nthe underlying data so that each new pixel is a weighted average of the\noriginal underlying pixels.  This greatly reduces the Moir\u00e9 patterns.\nHowever, when the *interpolation_stage* is set to 'data', it also introduces\nwhite regions to the image that are not in the original data, both in the\nalternating bands on the left hand side of the image, and in the boundary\nbetween the red and blue of the large circles in the middle of the image.\nThe interpolation at the 'rgba' stage has a different artifact, with the alternating\nbands coming out a shade of purple; even though purple is not in the original\ncolormap, it is what we perceive when a blue and red stripe are close to each\nother.\n\nThe default for the *interpolation* keyword argument is 'auto' which\nwill choose a Hanning filter if the image is being down-sampled or up-sampled\nby less than a factor of three.  The default *interpolation_stage* keyword\nargument is also 'auto', and for images that are down-sampled or\nup-sampled by less than a factor of three it defaults to 'rgba'\ninterpolation.\n\nAnti-aliasing filtering is needed, even when up-sampling. The following image\nup-samples 450 data pixels to 530 rendered pixels. You may note a grid of\nline-like artifacts which stem from the extra pixels that had to be made up.\nSince interpolation is 'nearest' they are the same as a neighboring line of\npixels and thus stretch the image locally so that it looks distorted.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fig, ax = plt.subplots(figsize=(6.8, 6.8))\nax.imshow(alarge, interpolation='nearest', cmap='grey')\nax.set_title(\"up-sampled by factor a 1.17, interpolation='nearest'\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Better anti-aliasing algorithms can reduce this effect:\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fig, ax = plt.subplots(figsize=(6.8, 6.8))\nax.imshow(alarge, interpolation='auto', cmap='grey')\nax.set_title(\"up-sampled by factor a 1.17, interpolation='auto'\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "Apart from the default 'hanning' anti-aliasing, `~.Axes.imshow` supports a\nnumber of different interpolation algorithms, which may work better or\nworse depending on the underlying data.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fig, axs = plt.subplots(1, 2, figsize=(7, 4), layout='constrained')\nfor ax, interp in zip(axs, ['hanning', 'lanczos']):\n    ax.imshow(alarge, interpolation=interp, cmap='gray')\n    ax.set_title(f\"interpolation='{interp}'\")"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "A final example shows the desirability of performing the anti-aliasing at the\nRGBA stage when using non-trivial interpolation kernels.  In the following,\nthe data in the outer circle is exactly 0.0, and data in the inner circle\nis exactly 2.0. If we perform the *interpolation_stage* in 'data' space and\nuse an anti-aliasing filter (first panel), then floating point imprecision\nmakes some of the data values just a bit less than zero or a bit more than\n2.0, and they get assigned the under- or over- colors. This can be avoided if\nyou do not use an anti-aliasing filter (*interpolation* set to\n'nearest'), however, that makes the part of the data susceptible to Moir\u00e9\npatterns much worse (second panel).  Therefore, we recommend the default\n*interpolation* of 'hanning'/'auto', and *interpolation_stage* of\n'rgba'/'auto' for most down-sampling situations (last panel).\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "a = alarge + 1\ncmap = plt.get_cmap('RdBu_r')\ncmap.set_under('yellow')\ncmap.set_over('limegreen')\n\nfig, axs = plt.subplots(1, 3, figsize=(7, 3), layout='constrained')\nfor ax, interp, space in zip(axs.flat,\n                             ['hanning', 'nearest', 'hanning', ],\n                             ['data', 'data', 'rgba']):\n    im = ax.imshow(a, interpolation=interp, interpolation_stage=space,\n                   cmap=cmap, vmin=0, vmax=2)\n    title = f\"interpolation='{interp}'\\nstage='{space}'\"\n    if ax == axs[2]:\n        title += '\\nDefault'\n    ax.set_title(title, fontsize='medium')\nfig.colorbar(im, ax=axs, extend='both', shrink=0.8)"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Up-sampling\n\nIf we up-sample, then we can represent a data pixel by many image or screen pixels.\nIn the following example, we greatly over-sample the small data matrix.\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "np.random.seed(19680801+9)\na = np.random.rand(4, 4)\n\nfig, axs = plt.subplots(1, 2, figsize=(6.5, 4), layout='compressed')\naxs[0].imshow(asmall, cmap='viridis')\naxs[0].set_title(\"interpolation='auto'\\nstage='auto'\")\naxs[1].imshow(asmall, cmap='viridis', interpolation=\"nearest\",\n              interpolation_stage=\"data\")\naxs[1].set_title(\"interpolation='nearest'\\nstage='data'\")\nplt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "The *interpolation* keyword argument can be used to smooth the pixels if desired.\nHowever, that almost always is better done in data space, rather than in RGBA space\nwhere the filters can cause colors that are not in the colormap to be the result of\nthe interpolation.  In the following example, note that when the interpolation is\n'rgba' there are red colors as interpolation artifacts.  Therefore, the default\n'auto' choice for *interpolation_stage* is set to be the same as 'data'\nwhen up-sampling is greater than a factor of three:\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fig, axs = plt.subplots(1, 2, figsize=(6.5, 4), layout='compressed')\nim = axs[0].imshow(a, cmap='viridis', interpolation='sinc', interpolation_stage='data')\naxs[0].set_title(\"interpolation='sinc'\\nstage='data'\\n(default for upsampling)\")\naxs[1].imshow(a, cmap='viridis', interpolation='sinc', interpolation_stage='rgba')\naxs[1].set_title(\"interpolation='sinc'\\nstage='rgba'\")\nfig.colorbar(im, ax=axs, shrink=0.7, extend='both')"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        "## Avoiding resampling\n\nIt is possible to avoid resampling data when making an image.  One method is\nto simply save to a vector backend (pdf, eps, svg) and use\n``interpolation='none'``.  Vector backends allow embedded images, however be\naware that some vector image viewers may smooth image pixels.\n\nThe second method is to exactly match the size of your axes to the size of\nyour data. The following figure is exactly 2 inches by 2 inches, and\nif the dpi is 200, then the 400x400 data is not resampled at all. If you download\nthis image and zoom in an image viewer you should see the individual stripes\non the left hand side (note that if you have a non hiDPI or \"retina\" screen, the html\nmay serve a 100x100 version of the image, which will be downsampled.)\n\n"
      ]
    },
    {
      "cell_type": "code",
      "execution_count": null,
      "metadata": {
        "collapsed": false
      },
      "outputs": [],
      "source": [
        "fig = plt.figure(figsize=(2, 2))\nax = fig.add_axes((0, 0, 1, 1))\nax.imshow(aa[:400, :400], cmap='RdBu_r', interpolation='nearest')\nplt.show()"
      ]
    },
    {
      "cell_type": "markdown",
      "metadata": {},
      "source": [
        ".. admonition:: References\n\n   The use of the following functions, methods, classes and modules is shown\n   in this example:\n\n   - `matplotlib.axes.Axes.imshow`\n\n"
      ]
    }
  ],
  "metadata": {
    "kernelspec": {
      "display_name": "Python 3",
      "language": "python",
      "name": "python3"
    },
    "language_info": {
      "codemirror_mode": {
        "name": "ipython",
        "version": 3
      },
      "file_extension": ".py",
      "mimetype": "text/x-python",
      "name": "python",
      "nbconvert_exporter": "python",
      "pygments_lexer": "ipython3",
      "version": "3.14.6"
    }
  },
  "nbformat": 4,
  "nbformat_minor": 0
}