.. DO NOT EDIT. .. THIS FILE WAS AUTOMATICALLY GENERATED BY SPHINX-GALLERY. .. TO MAKE CHANGES, EDIT THE SOURCE PYTHON FILE: .. "auto_examples/plot_0011_2d-image-interpolation.py" .. LINE NUMBERS ARE GIVEN BELOW. .. only:: html .. note:: :class: sphx-glr-download-link-note :ref:`Go to the end ` to download the full example code. .. rst-class:: sphx-glr-example-title .. _sphx_glr_auto_examples_plot_0011_2d-image-interpolation.py: 2d Image Interpolation =============================================================================== This example computes the optimal transport between two simple 2-dimensional images and then generates a simple approximation of the displacement interpolation .. GENERATED FROM PYTHON SOURCE LINES 7-43 .. code-block:: Python import time import matplotlib.pyplot as plt import numpy as np import scipy.sparse import MultiScaleOT # Create two images: a square and one rotated by 45deg # with approximately equal areas hierarchyDepth=7 # feel free to play with this value, up to 7 (i.e. 128x128 images) it should be quite low-weight n=2**hierarchyDepth nLayers=hierarchyDepth+1 # Square img1=np.zeros((n,n),dtype=np.double) thresh=int(n*0.5*(1.-1./np.sqrt(2))) img1[thresh:-thresh,thresh:-thresh]=1. img1=img1/np.sum(img1) # Diamond img2=np.abs(np.arange(n).reshape((-1,1))-n/2+0.5)+np.abs(np.arange(n).reshape((1,-1))-n/2+0.5) img2[...]=(img2zeroThresh) mu=img.ravel()[keep] pos=pos[keep] return (mu,pos) # extract measures from images mu1,pos1=extractMeasureFromImage(img1) mu2,pos2=extractMeasureFromImage(img2) .. GENERATED FROM PYTHON SOURCE LINES 62-63 Setup multi-scale solver .. GENERATED FROM PYTHON SOURCE LINES 63-82 .. code-block:: Python # generate multi-scale representations MultiScaleSetup1=MultiScaleOT.TMultiScaleSetup(pos1,mu1,hierarchyDepth,childMode=0) MultiScaleSetup2=MultiScaleOT.TMultiScaleSetup(pos2,mu2,hierarchyDepth,childMode=0) # generate a cost function object costFunction=MultiScaleOT.THierarchicalCostFunctionProvider_SquaredEuclidean( MultiScaleSetup1,MultiScaleSetup2) # eps scaling epsScalingHandler=MultiScaleOT.TEpsScalingHandler() epsScalingHandler.setupGeometricMultiLayerB(nLayers,1.,4.,2,2) # error goal errorGoal=1E-3 # sinkhorn solver object SinkhornSolver=MultiScaleOT.TSinkhornSolverStandard(epsScalingHandler, 0,hierarchyDepth,errorGoal, MultiScaleSetup1,MultiScaleSetup2,costFunction ) .. GENERATED FROM PYTHON SOURCE LINES 83-84 Solve .. GENERATED FROM PYTHON SOURCE LINES 84-91 .. code-block:: Python t1=time.time() SinkhornSolver.initialize() print(SinkhornSolver.solve()) t2=time.time() print("solving time: ",t2-t1) .. rst-class:: sphx-glr-script-out .. code-block:: none 0 solving time: 7.584433555603027 .. GENERATED FROM PYTHON SOURCE LINES 92-93 Extract coupling data in a suitable sparse data structure .. GENERATED FROM PYTHON SOURCE LINES 93-105 .. code-block:: Python couplingData=SinkhornSolver.getKernelPosData() # couplingData is a container for the coupling data in scipy.sparse.coo_matrix format # by calling the method couplingData.getDataTuple() one could obtain the list of # non-zero values and their row and column indices # we plug this into a simple routine for approximating the displacement interpolation at some time t t=0.5 interpData=MultiScaleOT.interpolateEuclidean(couplingData,pos1,pos2,t) # interpData is a container of particle masses and coordinates # these can be extracted via interpData.getDataTuple() muT,posT=interpData.getDataTuple() .. GENERATED FROM PYTHON SOURCE LINES 106-108 the intermediate measure can in principle be visualized as a weighted point cloud this can be slow on large images and also may not give a very good visual impression of the measure .. GENERATED FROM PYTHON SOURCE LINES 108-111 .. code-block:: Python plt.scatter(posT[:,0],posT[:,1],s=10*muT/np.max(muT)) plt.show() .. image-sg:: /auto_examples/images/sphx_glr_plot_0011_2d-image-interpolation_002.png :alt: plot 0011 2d image interpolation :srcset: /auto_examples/images/sphx_glr_plot_0011_2d-image-interpolation_002.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 112-117 alternatively, the point cloud can be rasterized to an image of suitable dimensions particle coordinates are assumed to be pixels (one pixel has length 1, origin at (0,0)) one may need to rescale coordinates appropriately before calling the function and one has to provide a suitable target image the target image is allocated: .. GENERATED FROM PYTHON SOURCE LINES 117-124 .. code-block:: Python reImg=np.zeros((n,n),dtype=np.double) # rasterize MultiScaleOT.projectInterpolation(interpData,reImg) # show rasterization plt.imshow(reImg) plt.show() .. image-sg:: /auto_examples/images/sphx_glr_plot_0011_2d-image-interpolation_003.png :alt: plot 0011 2d image interpolation :srcset: /auto_examples/images/sphx_glr_plot_0011_2d-image-interpolation_003.png :class: sphx-glr-single-img .. GENERATED FROM PYTHON SOURCE LINES 125-126 now do this for a whole sequence of times .. GENERATED FROM PYTHON SOURCE LINES 126-141 .. code-block:: Python nT=10 tList=np.linspace(0.,1.,num=nT) fig=plt.figure(figsize=(nT*2,2)) for i,t in enumerate(tList): fig.add_subplot(1,nT,i+1) # create displacement interpolations and rasterize them to image interpData=MultiScaleOT.interpolateEuclidean(couplingData,pos1,pos2,t) reImg=np.zeros((n,n),dtype=np.double) MultiScaleOT.projectInterpolation(interpData,reImg) plt.imshow(reImg) plt.axis("off") plt.tight_layout() plt.show() .. image-sg:: /auto_examples/images/sphx_glr_plot_0011_2d-image-interpolation_004.png :alt: plot 0011 2d image interpolation :srcset: /auto_examples/images/sphx_glr_plot_0011_2d-image-interpolation_004.png :class: sphx-glr-single-img .. rst-class:: sphx-glr-timing **Total running time of the script:** (0 minutes 11.849 seconds) .. _sphx_glr_download_auto_examples_plot_0011_2d-image-interpolation.py: .. only:: html .. container:: sphx-glr-footer sphx-glr-footer-example .. container:: sphx-glr-download sphx-glr-download-jupyter :download:`Download Jupyter notebook: plot_0011_2d-image-interpolation.ipynb ` .. container:: sphx-glr-download sphx-glr-download-python :download:`Download Python source code: plot_0011_2d-image-interpolation.py ` .. container:: sphx-glr-download sphx-glr-download-zip :download:`Download zipped: plot_0011_2d-image-interpolation.zip ` .. only:: html .. rst-class:: sphx-glr-signature `Gallery generated by Sphinx-Gallery `_