{"nbformat_minor": 0, "worksheets": [{"cells": [{"source": ["Sparse Spikes Deconvolution with Continuous Basis-Pursuit\n", "=========================================================\n", "\n*Important:* Please read the [installation page](http://gpeyre.github.io/numerical-tours/installation_matlab/) for details about how to install the toolboxes.\n", "$\\newcommand{\\dotp}[2]{\\langle #1, #2 \\rangle}$\n", "$\\newcommand{\\enscond}[2]{\\lbrace #1, #2 \\rbrace}$\n", "$\\newcommand{\\pd}[2]{ \\frac{ \\partial #1}{\\partial #2} }$\n", "$\\newcommand{\\umin}[1]{\\underset{#1}{\\min}\\;}$\n", "$\\newcommand{\\umax}[1]{\\underset{#1}{\\max}\\;}$\n", "$\\newcommand{\\umin}[1]{\\underset{#1}{\\min}\\;}$\n", "$\\newcommand{\\uargmin}[1]{\\underset{#1}{argmin}\\;}$\n", "$\\newcommand{\\norm}[1]{\\|#1\\|}$\n", "$\\newcommand{\\abs}[1]{\\left|#1\\right|}$\n", "$\\newcommand{\\choice}[1]{ \\left\\{  \\begin{array}{l} #1 \\end{array} \\right. }$\n", "$\\newcommand{\\pa}[1]{\\left(#1\\right)}$\n", "$\\newcommand{\\diag}[1]{{diag}\\left( #1 \\right)}$\n", "$\\newcommand{\\qandq}{\\quad\\text{and}\\quad}$\n", "$\\newcommand{\\qwhereq}{\\quad\\text{where}\\quad}$\n", "$\\newcommand{\\qifq}{ \\quad \\text{if} \\quad }$\n", "$\\newcommand{\\qarrq}{ \\quad \\Longrightarrow \\quad }$\n", "$\\newcommand{\\ZZ}{\\mathbb{Z}}$\n", "$\\newcommand{\\CC}{\\mathbb{C}}$\n", "$\\newcommand{\\RR}{\\mathbb{R}}$\n", "$\\newcommand{\\EE}{\\mathbb{E}}$\n", "$\\newcommand{\\Zz}{\\mathcal{Z}}$\n", "$\\newcommand{\\Ww}{\\mathcal{W}}$\n", "$\\newcommand{\\Vv}{\\mathcal{V}}$\n", "$\\newcommand{\\Nn}{\\mathcal{N}}$\n", "$\\newcommand{\\NN}{\\mathcal{N}}$\n", "$\\newcommand{\\Hh}{\\mathcal{H}}$\n", "$\\newcommand{\\Bb}{\\mathcal{B}}$\n", "$\\newcommand{\\Ee}{\\mathcal{E}}$\n", "$\\newcommand{\\Cc}{\\mathcal{C}}$\n", "$\\newcommand{\\Gg}{\\mathcal{G}}$\n", "$\\newcommand{\\Ss}{\\mathcal{S}}$\n", "$\\newcommand{\\Pp}{\\mathcal{P}}$\n", "$\\newcommand{\\Ff}{\\mathcal{F}}$\n", "$\\newcommand{\\Xx}{\\mathcal{X}}$\n", "$\\newcommand{\\Mm}{\\mathcal{M}}$\n", "$\\newcommand{\\Ii}{\\mathcal{I}}$\n", "$\\newcommand{\\Dd}{\\mathcal{D}}$\n", "$\\newcommand{\\Ll}{\\mathcal{L}}$\n", "$\\newcommand{\\Tt}{\\mathcal{T}}$\n", "$\\newcommand{\\si}{\\sigma}$\n", "$\\newcommand{\\al}{\\alpha}$\n", "$\\newcommand{\\la}{\\lambda}$\n", "$\\newcommand{\\ga}{\\gamma}$\n", "$\\newcommand{\\Ga}{\\Gamma}$\n", "$\\newcommand{\\La}{\\Lambda}$\n", "$\\newcommand{\\si}{\\sigma}$\n", "$\\newcommand{\\Si}{\\Sigma}$\n", "$\\newcommand{\\be}{\\beta}$\n", "$\\newcommand{\\de}{\\delta}$\n", "$\\newcommand{\\De}{\\Delta}$\n", "$\\newcommand{\\phi}{\\varphi}$\n", "$\\newcommand{\\th}{\\theta}$\n", "$\\newcommand{\\om}{\\omega}$\n", "$\\newcommand{\\Om}{\\Omega}$\n"], "metadata": {}, "cell_type": "markdown"}, {"source": ["This numerical tour explores the use of the\n", "continuous Basis-Pursuit (C-BP) method to perform sparse deconvolution.\n", "This method was proposed in [Ekanadham11](#biblio)."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 2, "cell_type": "code", "language": "python", "metadata": {}, "input": ["addpath('toolbox_signal')\n", "addpath('toolbox_general')\n", "addpath('solutions/sparsity_9_sparsespikes_cbp')"]}, {"source": ["Sparse Spikes Deconvolution\n", "---------------------------\n", "We consider the problem of estimating an unknown Radon measure $m_0 \\in \\Mm(\\mathbb{T})$ from low-resolution noisy observations\n", "$$\n", "      y=\\Phi(m_0)+w_0 \\in L^2(\\mathbb{T})\n", "$$\n", "where $ w_0 \\in L^2(\\mathbb{T})$ is some measurement noise, and\n", "$ \\Phi : \\Mm(\\mathbb{T}) \\rightarrow L^2(\\mathbb{T}) $ is a convolution operator\n", "operator with a smooth kernel $ \\phi \\in C^2(\\mathbb{T}) $, i.e.\n", "$$\n", "      \\forall x \\in \\mathbb{T}, \\quad (\\Phi m)(x) = \\int_{\\mathbb{T}} \\phi(x - y) d m(y).\n", "$$\n", "We focus our attention here for simplicity on the compact 1-D domain\n", "$\\mathbb{T}=\\RR/\\ZZ$ (i.e. an interval with periodic boundary conditions)\n", "but the continuous basis-pursuit method can be extended to higher dimensional settings.\n", "\n", "\n", "In order to recover sparse measures (i.e. sums of Diracs),\n", "it makes sense to consider the following regularization (sometimes called\n", "BLASSO for Beurling LASSO in [deCastroGamboa12](#biblio))\n", "$$ \\umin{m \\in \\Mm(\\mathbb{T})} \\frac{1}{2}\\norm{y-\\Phi(m)}^2 + \\la \\abs{m}(\\mathbb{T}) $$\n", "where $\\abs{m}(\\mathbb{T})$ is the total variation of the measure $m$, defined as\n", "$$ \\abs{m}(\\mathbb{T}) =\n", "  \\sup \\enscond{ \\int_{\\mathbb{T}} g(x) d m(x) }{ g \\in \\Cc(\\mathbb{T}), \\: \\normi{g} \\leq 1 }.$$\n", "This formulation of the recovery of sparse Radon measures has recently\n", "received lots of attention in the literature, see for instance the works\n", "of [CandesFernandezGranda13,deCastroGamboa12,DuvalPeyre13](#biblio).\n", "\n", "\n", "The BLASSO optimization problem is convex but infinite\n", "dimensional, and while there exists solvers when $\\Phi$\n", "is measuring a finite number of Fourier frequency\n", "(see [CandesFernandezGranda13](#biblio), and the related Numerical Tour\n", "\"Sparse Spikes Deconvolution over the Space of Measures\"),\n", "they are do not scale well with the number of frequencies.\n", "The vast majority of practitioners thus approximate the recovery over arbitrary\n", "measures by a finite dimensional problem computed over a finite\n", "grid $z = (z_i)_i \\in \\mathbb{T}^N$,\n", "by restricting their attention to measures of the form\n", "$$ m_{a,z} = \\sum_{i=1}^N a_i \\de_{z_i} \\in \\Mm(\\mathbb{T}).\n", "$$\n", " For such a discrete measure, one has $ \\abs{m_{a,z}}(\\mathbb{T}) = \\norm{a}_1$,\n", "which can be interpreted as the fact that $\\abs{\\cdot}(\\mathbb{T})$ is the natural extension of\n", "the $\\ell^1$ norm from finite dimensional vectors to infinite dimensional space of measures.\n", "Inserting this parameterization inside the BLASSO leads to the celebrated Basis-Pursuit problem\n", "[Chen99](#biblio),\n", "$$ \\umin{a \\in \\RR^N} \\frac{1}{2}\\norm{y-\\Phi_z a}^2 + \\la \\norm{a}_1.\n", "$$\n", "where in the following we make use of the notations\n", "$$ \\Phi_z a = \\sum_{i=1}^N a_i \\phi(\\cdot-z_i),\n", "$$\n", "$$ \\Phi_z' b = \\sum_{i=1}^N b_i \\phi'(\\cdot-z_i).\n", "$$\n", "One can understand the BP problem as performing a nearest neighbor interpolation of the Dirac's locations.\n", "\n", "Continuous Basis-Pursuit\n", "------------------------\n", "To obtain a better approximation of the infinite dimensional problem, [Ekanadham11](#biblio)\n", "proposes to perform a first order approximation of the kernel.\n", "This method assumes that the measures are positive.\n", "\n", "\n", "To ease the exposition, we consider an uniform grid $z_i=i/N$ of $N$ points, so that\n", "the grid size is $\\De=1/N$.\n", "The C-BP method solves\n", "$$ \\umin{(a,b) \\in \\RR^N \\times \\RR^N} \\norm{ y - \\Phi_z a - \\frac{\\De}{2} \\Phi'_z b }^2 + \\la \\norm{a}_1 $$\n", "subject to\n", "$$ \\abs{b} \\leq a \\qandq a \\geq 0\n", "$$\n", "where the constraint inequalities should be understood to hold\n", "coordinate-wise.\n", "This is a convex optimization problem, which can be solved using traditional\n", "conic optimization methods. We detail in this tour\n", "a simple proximal splitting scheme which can be used.\n", "\n", "\n", "If $(a^\\star,b^\\star)$ are solutions of C-BP, one recovers\n", "an output discrete measure defined by\n", "$$ m^\\star = \\sum_{a^\\star_i \\neq 0} a^\\star_i \\de_{x_i^\\star} \\qwhereq x_i^\\star = z_i + \\frac{\\De b_i^\\star}{2 a_i^\\star} .\n", "$$\n", "The rationale behind C-BP is to perform a first order Taylor approximation of the\n", "kernel $\\Phi$, where the variable $\\tau_i =  \\De b_i /(2 a_i) \\in [-\\De/2,\\De/2]$\n", "encodes the horizontal shift of the Dirac location with respect to the grid sample $z_i$.\n", "The landmark idea introduced in [Ekanadham11](#biblio) is that, while the optimization\n", "is non-convex with respect to $\\tau$, it is convex with respect to $b$.\n", "\n", "Convolution Operators \n", "----------------------\n", "\n", "\n", "\n", "We define continuous kernel, defined on $ [0,1] $, and select here\n", "Gaussian kernels\n", "$$ \\phi(t) = \\exp\\pa{ - \\frac{t^2}{2 \\sigma^2} },  $$\n", "$$ \\phi'(t) = -\\frac{t}{\\sigma^2} \\exp\\pa{ - \\frac{t^2}{2 \\sigma^2} }.  $$"], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 3, "cell_type": "code", "language": "python", "metadata": {}, "input": ["sigma = .03;\n", "phi  = @(t)exp(-t.^2/(2*sigma^2));\n", "phi1 = @(t)-t/(sigma^2).*exp(-t.^2/(2*sigma^2));"]}, {"source": ["Set the discretization size $\\De = 1/N$."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 4, "cell_type": "code", "language": "python", "metadata": {}, "input": ["N = 32; \n", "Delta = 1/N;"]}, {"source": ["We discribed the method using continuous observation $y \\in\n", "L^2(\\mathbb{T})$. In fact, it can be applied to a discretized obervation\n", "space as well. We set $P$ to be the number of observation points."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 5, "cell_type": "code", "language": "python", "metadata": {}, "input": ["rho = 64;\n", "P = N*rho;"]}, {"source": ["Useful helpers: computer periodic convolutions, and adjoint of convolution\n", "(i.e. convolution with the time-reversed kernel)."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 6, "cell_type": "code", "language": "python", "metadata": {}, "input": ["convol  = @(x,h)real(ifft(fft(x).*fft(h)));\n", "convolS = @(x,h)real(ifft(fft(x).*conj(fft(h))));"]}, {"source": ["Up-sampling and down-sampling operator from the sampling to the observation\n", "grid."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 7, "cell_type": "code", "language": "python", "metadata": {}, "input": ["upsample = @(x,rho)upsampling(x,1,rho);\n", "downsample = @(y,rho)y(1:rho:end);"]}, {"source": ["Observation grid and sampling grid."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 8, "cell_type": "code", "language": "python", "metadata": {}, "input": ["t = [0:P/2, -P/2+1:-1]' / P;\n", "t1 = (0:P-1)'/P;"]}, {"source": ["Discretize the filter on the observation grid."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 9, "cell_type": "code", "language": "python", "metadata": {}, "input": ["phi_d = phi(t);\n", "phi1_d = phi1(t);"]}, {"source": ["Operator $\\Phi_z$, $\\Phi_z^*$, $\\De/2 \\Phi_z'$, $(\\De/2 \\Phi_z')^{*}$,\n", "$ \\Ga_z = [\\Phi_z, \\De/2 \\Phi_z'] $ and $\\Ga_z^*$."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 10, "cell_type": "code", "language": "python", "metadata": {}, "input": ["Phi  = @(x)convol(upsample(x,rho),phi_d);\n", "PhiS = @(x)downsample(convolS(x,phi_d),rho);\n", "Psi  = @(s)convol(upsample(s,rho),phi1_d)*Delta/2;\n", "PsiS = @(s)downsample(convolS(s,phi1_d),rho)*Delta/2;\n", "Gamma  = @(u)Phi(u(:,1)) - Psi(u(:,2));\n", "GammaS = @(y)[PhiS(y), -PsiS(y)];"]}, {"source": ["Generate the input signal $m_{a_0,x_0}$ composed of three spikes.\n", "The variable $\\kappa \\in [0,1]$ controls the amount of shift between the locations\n", "$x_0$ and the sampling grid. The larger $\\kappa$, the worse is the\n", "approximation performed by C-BP. Setting $\\kappa$ close to 1 will make\n", "C-BP fail to estimate correctly the number of spikes of the input\n", "measure."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 11, "cell_type": "code", "language": "python", "metadata": {}, "input": ["k = 3; % number of spikes\n", "kappa = .9;\n", "I = round( N/(2*k):N/k:N );\n", "a0 = zeros(N,1); a0(I) = [.6 1 .8];\n", "d0 = zeros(N,1); d0(I) = [-.2 1 -.7] * kappa;\n", "b0 = d0.*a0;\n", "x0 = (0:N-1)'/N + d0*Delta/2;"]}, {"source": ["Generate the true observations $y=\\Phi(m_{a_0,x_0})$. Note that here we\n", "assume there is no noise, i.e. $w_0=0$. It is important to realize\n", "that, even if there is no noise, when the spikes are off the grid, the\n", "linear approximation made by C-BP (and of course also by BP) create an\n", "error term which should be interpreted as noise. This thus necessitate\n", "the use of a well-chosen value of $\\lambda$, which should scale like\n", "$\\De^2$ (magnitude of the approximation error)."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 12, "cell_type": "code", "language": "python", "metadata": {}, "input": ["y = zeros(P,1);\n", "for i=1:length(x0)\n", "    T = t-x0(i); T = mod(T,1); T(T>.5) = T(T>.5)-1;\n", "    y = y + a0(i) * phi( T );\n", "end"]}, {"source": ["Generate the basis pursuit approximation (quantification) $y=\\Phi_z\n", "a_0$."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 13, "cell_type": "code", "language": "python", "metadata": {}, "input": ["y0 = Phi(a0);"]}, {"source": ["Generate C-BP approximation (linear approximation) $y=\\Ga_z\n", "(a_0,b_0)$."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 14, "cell_type": "code", "language": "python", "metadata": {}, "input": ["y1 = Gamma([a0 b0]);"]}, {"source": ["Display."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAAcC0lEQVR4nO3d7ZKj\nONYuUDgx9/3OXDnnh7toCjDmU9pbrBUdHVlOZ1qJQI8kQPTDMHQAkM3/q10AADhDgAGQkgADICUB\nBkBKAgyAlAQYACkJMABSEmAApCTAAEhJgAGQkgADICUBBkBKAgyAlAQYACkJMABSEmAApCTAAEhJ\ngAGQkgADICUBBkBKAgyAlAQY7NX/Mf7z0c8q8CmQ2n9qFwBy6Pt+GIbl108r9kGQTrnjEPJaJtZs\nYPT57vTF2SvT8Fu+MgzD6kjr8/r0u8vfM75t9tHQPCMwOG85JpuGxydUxiQbvzV9ZfmDs9e774O/\nbxkJLyHA4DY/h0obLp7rWv1oaJsAg99Wh0TLyJm9p9sdY8sfPKrWKTqoSIDBLtujq3Ea8NupqWWi\nLENx9obp7OLP0dW3s2jQMD01AFJyHxgAKQkweJzJPXiCKUR4yrdbu8qXBJrU4AhMb5cIvu2H9k+4\nS1NXIWoaAN6jqRHYMAzmZwhp6Lr5Ch3ARU2NwJa0FNQ2/P21NeZJI/54oPEA62LUQZCVESIUI0IZ\nChTjTz4tP+KfDFsukFhLhDIoRsBipOhjNTWFCEEMw/B3evXj2Gst1YAzBBg8bdmTlWFwgwYDrPrQ\nGyazL9P0SjAnA4k0GGAQ1ti5ynB+AaKrf6rwURHOhfI2Yzit7nrb34UgUjSeRmBQVPg2AdIQYHCn\n/QMss4hwkQADICUBBqWZRYRbCDC4zbf5w77vV9c1MIsIVwgwAFJqfy1ECGgYDL8iGgfK8a8gpxNg\ncJdzN3j1vVNiIczmeD//FGPBmUIE3s7js5MSYFCHzj1cJMDgTidiSS+/rr+HWYPHZyfiHBjcYLuV\ncyolg2Htn9IrNCMwqEauhfGtJtRQaAIMbiOQMvp7fNz//fhso+fQBBhcdf0sifMsFW0/fVTVRCbA\ngPf6kl5/vSLDwhJg8LhvayF2Zh3DGIZhOls4+ycxCTC4x8XmTje/vOXiKbPQGv+ldmISYHCJpg1q\nEWDAq20PnQ3CIhNgUJlTLVUIpAYIMLjBLSGkSYVDBBicJ3KSOvTsG7OIYVkLER7nguxEVFYiAqym\n2b1BjpzX8oDmKvYfcCooJlOI1SzvbPXghqRu7HjYBQq4spFVUCgCrA5PgG2AuoK6BBjwRkfHzSb4\nAxJgFUyGWcPf/82+SyM21kIcaR/L+Hl4/awsB2gcAqyiZYulDcvn9uDRPsJOAqwWT4DNTcykdq7b\nYZQcjQCrYusJsDIMHnJXt0P3JQgBVtrfD9BbeQJs5/B4Kx18OESAFfUrmWRYMg9Fjtp/1JVa08kI\nxUocdQxDN7vscFyGQ+MVnzrK6N5a63thVl++ABuvcF0uvDS9+DXgskzf1g9dPgT2805HSDMC7o18\no7ISSRZgfd9PRir9clez85GaNfeedr2FUEdxtHYObM8do1V4fEN7Hu0sqfp7PbE91VF1yUZgP31G\nYLOB2vINcJpmi1bF7P1vaCrAVsMpQmIdGn6N73QmDO5116HU6izitLVMEWbtTCGm2Nyrwk57cpdD\nVay/crtDh9ehynLg1pVsBDYMw/IqxM+E4eq3Ijgx/BrfbxAWVoEaUe+wLVmAdWvhNL4SKrdokh43\nxNHOFGJMp4dfV34KGF08Br9xqXAEAiwHB8mr6LjAHgKsBO1RY4pVqI4LbBBgD7ql9TFTEYcqyOuJ\nPoeOaXX5LuJoj2tPmneiilu906iwE9vwRGW5XrQWI7Cn3Hjq2CAsFE0VBCHAYJdavQe9luue63Po\nzdQlwJ5149I1wCEls18/owoB9ojn9mbHyXvotcA2AVaftRDju3hG82IV2zuuOFplRytLP6MiAXY/\nd/5DdeUPEwdmeQIM4tK7hw0C7GYPDb+e+53sV3H7690f8uhhuGR2pBYBlpLjpCRbG2ISYI8wVOIu\n9iX4RoDdqUBX3WRFLdWDRI3vVHj+cPZZqqkkayHWZy3EyG5akVkVp6GyEjECu02xfp++3ttoUU8o\nv9FUU3kCDL6qMhm1QZflpwibKEIZXkKA3cNluwCFCTBIIMgQMLjqI2Y9y8IE2A2qLJTnUHnaja3h\njctdqu6nWZs0EQEGtKD68Gv26UKwAAF2VcXDxqHynCCt4ZTqhhkBBqQXqsOhq1GMALuk+mHjUHlC\n9Wr9RnXDlAA7TyMCEQTscOhqlCHAbuCkcUsCtoZTqhtG1kI86cZm7vraa8PwT3n6Pmiz+3KW13vO\n7R2OuyrLUVmAEdgZkTu/kcsWX/Dh14dB2CjLFshSznQE2CVxmrk4JckrUSsjw2Zi7v8xS9USAXZY\n2E66Ru0u0Wp222urO+yROOWofJQAOyb4MeNoOS14zS5lKedDMu7hGcscnAA74KH976G11xwt+z2d\nXg9V8Wv7K9O/9/Yqu72ypiV8W009TYDt0vfPHjM3mh0tDpifUm+iaYal/kP2y3IkTsmwhwiw32YH\nzIljpvDi1rNCjh8eYY3tCGXo/hSjer/klq1xsXEMVSO/3vN4fT23NQ71LINUSnxN3Qc21vpdd3LM\n9qIs3b2P8TaUTqdvxZCxI//Nsq6z/0UzqY/E0bSaukZrqrB2Aqzv+zG3pl8f/z3rr2fcz2YHzLTV\nzvjnXLdauW1sitXGsUv+17VXX5/Ct1dTtbQTYPt9bozfPyhJvWMtD5iPjZxuY+GAPfXbwJ85tVrX\nwSt6cjAOP6usemnvsqOm/uprRqipmM6PVKJZHYFdmzprct6tkeq+psmaXWqmrtuur7jVFD8c3jgC\nW/pST+Fr7w5vOD22qN9X1OxMoop+eX0lqqnqGg+w+D2I6myil1DRWaip/dqZQuweuAoRgLCaCjAA\n3sONzACk1PI5sIozitsffeU2tVvKUHLLfPuswrWTpUZWv1umGBFqZLYCRa1NMf3ue46RjJoNsLvu\na773o4utEPPzz/9zm8GzW2a7GGXK8LMYZSplZ43UKsa0sa5YIyUPlp97Rd1jpGILlogpxKKGYYiw\nIyrDKE7T8NCK9SfKEGGDBCkGwTU7AuOn6m1E9fY6lGKjn+0C1C1DKMWmEPcUg1UC7I2CHJljk12r\nAH/WpP/n/xU3SPW6CKV6gkaYvhuGYTxOZdg3phBfqnoDUfHTR8MfXdUNEmRrEMonOPVstrU8VxDn\nKsTlpRxVLm36fG71C71mxahyzVuoGln9VvliBKmRkiMeNZJdywEGQMNMIQKQkgADICUBBkBKAgyA\nlAQYACkJMABSEmAApCTAAEhJgAGQkgADICUBBkBKAgyAlAQYXNX/MX1l9sW5X3u1ZNA0D7SESyI8\n/BDeSYDBzWaPeF59qtPyaVjdl0dUrz4sqvOMKBBgcLvPowjHgFnG0uqgbfk4x3Facvkt6QWdAIP4\nnAyDVQIMQpvNHAIjAQaXfCYMx6/H1zcu6Pj2I7Nvbb8TcNEUACm5DwyAlAQYACkJMABSEmAApCTA\nAEhJgAGQkvvAoKbp3cluaYFDGg8wixcQ21+R9VkEuFJJYC7+XcKNB1gXow6CPGUjQjEilCFIMb50\nrkqXK8KmUIyAxUjR+3cODCqYzRx+/lt+C9ggwKCm1X62DIM9BBiUNn002PT12pNGkEz9mdZHRZhK\nhpkxwLaHX/bcWlKc/nnCzscpxNH+RRwQU/jG4dXit923yxjbphChqJ+thKs5YCcBBkBK+QKs/+PQ\ntyCU7Qmq901f0XWLFmz5xbnfeUPJokp2Dmx6XnH1HOPnlRSnH3mho41J3wuzt/jZuLGULMC2qXKg\nDZ+R05hk40BqdqHg7MXPK8MwLAdeszeP73zsLyihqQD7mHVeZhWZvcJ4iWH4Z7hmEBbBc/Nw3yr3\nE0JjezWLse77iG32g2MQdosJqmVjmG6+sakAW+1TSCyisUtSxc982uj6x9RUgHXiisCONggGYXE0\nsP1nY7I2JAuw6dzubPg8HSl3kgxIZaNx2//+5bd+vjm1xq91cTEPcXwakKP747mf4op3thvLqwfi\nb4R894HBCzU06wO3STaFCElJoAgslNwYIzAIzdKIt+j7vzbg7J8kJcCAlm1klQzLToBBOedmrgzC\nTptusWH497/VN5COAIPHaSWrmKXXVORzYN/u0zp3/1ZLd30tCTBIIHKDG9/q1gs4rvUwjaNchQiZ\nWJVjp0MXHP7cqv1/n8qV4b//fvB0wcNpkk0XathYzHe24OFf5f++7O/df1BR+QLs5/3kKe6/44Xs\nlcXsTK9xsa6Ypu3YbIne7stivsvnSY2juuW3Gmgnk7X124/MObQQCxRz12oabmPaY+PU1/b7J5OK\nNduNb49Q+Rlg0x/cWPnw20MTM67EkW8EtmF17Ax12R9rCd/8brn3+ZZNruTbNRZgqzwPjGZYn/6n\nZgap3x5KOXvu1/j1xs/uX8k3Xby1H2ASC17iXHpF6xZ8O0e1bMo23nD0Z5cvpggzl9FDCXe1jAEv\n/g7CBnmhZAH2GQt/tDqrC1xxoq+gW5BUvinEQ+NoqOuJBjHalFcEzZz64pB8AQbIsKlb0mt2T5h5\nnRQEGJDY7UHT92Zx0kh2DgwyeqI9dNqmO37PMo0xAoOu+77Y3XSpumJ2FublE4m3p1fwlaVYSrBY\nyBUpVkOhop2LtJ6LsRPnZvaUZ1qY1w5BHvrDXQwyStF4JijiFeXrYKMH1/SWzmc1Kv7Khr/fcCLD\njraGG5+4UdoXZthzf7IAGwmw+srUwdFph6Y3eQLLMNgIp+mbj2bYsSd67PigbzH2qgx7+o+9a+Xl\n7ARYfY/Wwc/c+nzyt7c1veGDOhRdqz91KMP2B9jRj1hJsv/NT4+1Z3Yorf6N/X/7i6ctbxyEXS9M\nRQKstAKPU9kOrd/t1OLHG9r8cf2cLdz/G04E2O3ptfqD/5okWTN712p0HX3I5M5teyXADhUpeLYJ\nsKJWnz5wvQ52DrOu/M5WaiCWjabk5BUZxzPsaIBdLNiKP2GWcR+bH3r/d88FgrsGuMdnEW95ZHOo\nSBNgRe0MsNmuuX6Q/G9tm/xfv6dv+9mP/zm7vjmBcDTGPr9t/P3TD9rzg4e/FfLK7GmBH+qAb3z0\noV+1syN/Mb1Wf9Ve/xu6z7Xj/+27/5VrCcZj8PO5/7w6HnqngmrnWcw97+921925338x6v7a/5+s\nNQFW1HqA3dEtIrUbe7WHwmZPI3jlCpGtj259t79lqLqnc3lo+veugsURakS4SoDRitVxc9d1994n\ntDvDDgVYgZYi0bHwmWB4aJvsjJyfs4i3RNfqh/7rplnT0wRYOQ+dAyvjrrNiOy/oz7BJTnr6Prz9\nY6b9LWD8ZqIxe7Jnu/9xb3rFvHk0RePZzlJS+x+bHdDsgvvDN8BeuzayJdM/drZZbjmrN559JK/x\nFPXnn/tPJ3dFoutVB+xFCTL2ihSdiJmdl9oLrf1uv/X15+DpZxfE8CuC7TRajqEfTa9ox2yKxrOd\nEVgzlrc/mxi8aLpJI6yzYAwXxJ6h2GfgfvsZr1ctnvKcBBl7RYpOxIY90ZX57yvt3kUW/vlVF86g\nGH4FUfgG5BTplaLxNAILLfz+k0zh549IryxmQ7Htt12UIr2ySJCxV6ToRFDYXeOwjSjamKgUYPHN\nk+y++4UTrXafovFMUMQrUtQB5T2aYRu/XHrlcm/e5Bp7pWg8/1/tAkBNnsDLhjtvgU+VXlkIMN7o\nrhbk0MKMhl95XezoSK+HCDBe6uuCzmdtZ5j0ei3p9Zx8Adb/sfGGkuXh5b4FkqaqDTd2dOwSt0sW\nYJ/zih/LoNoONpi5q236OZFo+PVaiS47zChZgG37BFvtUvBqywyTXq+lO/20BBdKTq0uOb/9ntl3\nc/29FHD7JfX//La/l/2VXqmd20nSDb9mDWb81jJ0gC235okAi/wHEsET60vNSK/sTuwk6dJrJkXj\nmaCIUwKMJ9yZYf38OYTSqw2HdpLs6dUlaTwTFHFm+dCv2YYWYBx1c4Blbrb4Zv9O0sZ18ykaz3wX\ncYxXIU5fmb2heKHI7a7LEZ20b9jOnaSN9MoiX4ABxCe9ChBg0HUPLMxBe37uJA2c+spFgMHNNF5v\nsMww6VWeAIN/XByEGbq9wTSc+v7f56NKryoEGMABs4hy1UZFAgz+pQFij2FY2VXsPOX9p3YBIKK+\nP9keacXeQ11XZwQGQEr5RmDLlTj2fAt2GoZ/z8yfWPgOKCZZgP1cJurzSopFUAC4oqkpRKHFLU7v\nR3ZAKCnZCGyP5dq+0+8KOfbbOYto/pA2pHuifegAO5o9n/db2xfghI2uf0yhA+xE9ogrbnHuUg6g\npNABtvR5KPP49eeLz5zh53UXIlKLPQ4KSxZg3VoyfV6RWNxr5yAsw0QLtKmpqxABeA8BBl/tX5/e\n+B/KE2BwnvlDqEiAlRDkgtQIxYhQhu5IMR59UnOErRGhDJ1i/C1IMeITYHCSZxhCXQIMfnh0EAac\nJsDggDHDDL+gOgEGv01Tqu+lF4TQ+GNHnAvlVsuDxQ5Gs+KnQ+MBBkCrTCECkJIAAyAlAQZASvlW\no9+v4qNVtj969szo8mUouWW+fVbh2slSI6vfLVOMCDVS+PnpP/eKAmXYKIaHQ/3UbIBNm6QyzdOe\njy52VeTPP//zytNbZrsYZcrwsxhlKmVnjdQqxrSxrlgjJQ+Wn3tF3WOkYguWiCnEooZhiLAjKsMo\nTtPQ9331uz4+ZYiwQYIUg+CaHYHxU/U2onp7HUqx0c92AeqWIZRiU4h7isEqAfZGQY7MscmuVYDP\nR4//r7hBqtdFKNUTNML03TAM43Eqw74xhfhS1RuIip8+Gv7oqm6QIFuDUD7BqWezreW5gjhXIS4v\n5ahyadPnc6tf6DUrRpVr3kLVyOq3yhcjSI2UHPGokexaDjAAGmYKEYCUBBgAKQkwAFISYACkJMAA\nSEmAAZCSAAMgJQEGQEoCDICUBBgAKQkwAFISYACkJMAASEmAwT02Hut14olfHhIGP3kiM1wlbKAK\nAQZXzR6BOM2z8bGEy0doLn929eHxq0877DzkEAQYPGGaLp9Y+vYA6PHr2dvG2OvWniItvaATYHC7\n2x8Gb4oSVgkwuNnqAOu02ZgMGAkwuNm3U1nTCcDVIdryB28fzEFLrnYPAaAK94EBkJIAAyAlAQZA\nSgIMgJQEGAApCTAAUmr8PjD3fnLA8o6SirtPqMLwSvFvsmo8wLoYdXB9OYZmihGhDKvF6FfzYeiG\nlSR5sBjlCxO2RhSjehnqFmAPU4jwV2AM3fBoaB0yK8x6sMFbCTD415gW4xflM2P8xGVhgCkBxtst\nA6Oib4WpGKgQlgDj1TbSK2xmRCsP1CLASqh+PvYjQjEilKHLWYyHBogZN8VzFCMXAcZ7/Zw8LDyp\nuL88BmHQCTBe69CpL4EBAQkwiGU7UA3CYCTAeKP9w69is4gnAkmG8XICDHaJkxYRLveHCAQYrxPq\nxq9zTCRCJ8B4m+Atft5MhfIEGC+1PyoKhMqJWDUIAwHGi1ycPAwbFWELBo9q/3Eq0KqhG0TX1HJr\nVJySnRXG5PATjMB4i/jXbpwomInEj77rV7dAlc2yWpjPiy+vptsJMF7hYno9+1jLmxq11zaO2394\n4dj4+VmvraYnmEKkfTc2GX3XRxvATScSAxbvUT+n6QqnxezJqN+++/niVTX1ECMwGrfdprThnU9t\n3vMc7VpTrNuF+XhPTT1HgNGyROl1sXivyrDZrGCQ1SP3TFPPgtZZsYsEGM26N72C51/3jgxbRteh\nenlusxw6yWoodhcBRpuyjL3ubbzazrDTF6YH3AEMxW4hwGjNofmlc7//3l94ryabxYsDr+7hicTT\n17iKsYuaugqx7//sRh7I/Uqh7mOtaHaDc9Jr3r415en+kJ9W66tr8S+9XTsB1vf9mFvTr+/8iK4f\nd7XZPvf5Z8Un0C93+jdcUb3RXW3+b9/2+fNXm8XpG0LZHnzcchbz8xH3Hhq35M2yvpb/DFhl1bUT\nYOuGm2cMZndyzF5/6Ag8dGvkRpv1EhmP84fKvNosfiTaNzJW6Dnb85ylqyzDVm89wCJZ3f8sBXSL\nMm1c0rUHM17Z8Vyi3zsIe2i6L2OVVdF6gPWVz4fdvq7Me3qjb1C+bbL/pFOryvq+j7+ztB5gta3u\nfDubLW1NWG84v9iqG4fRrraorp0AG4Yhy1WIdneoTi+kAe0EWBc+t2CVZjQ11VeRG5nhAK1VA265\ncsq1FREIMABSEmBQhy58RRdH0i7fCEKAwRnipw3qMTUBBnCG4Vd1Agxq0gjWcvpSDoO2OAQYHCNy\nIAgBBhXoxUdwoi/i8o1QBBjwdvoTSQkwOEmr91qGX0EIMOC9Dl3KocsSjQCDanTk4QoBBocJnpbs\nHIS5fCMgAQZASgIMSnMqJZqfgzDDr5gEGJwnit5ALYclwAB2nQkz/IpGgAF03ZcMM3kY2X9qFwBe\nSoMYmWnDFIzA4Azx06Rv1aq6Y0owAuv7P0P4Ydh4ffzn8p0Qh659cJ+sMnOYQvQA6/t+TKOfX8st\nyuu7XhvXHnWaQlNTiH3fT8dhADQs+gjskHEucToUm0WaURrAqnQDgFgBdiVsvr1ZYhGQGSoC2uj6\nxxQrwE6HzWzUBQUM3eCKDKgoVoAtDcOwvNpwGIbV16f/hIAEHtwoeoB1i0Aa//ntdQDeoKmrEKEK\n4yqoQoABkJIAg9Jcggi3EGBwniiCigQYFOJUGdxLgMENhBOUJ8AASEmAAZCSAIOiXPcBdxFgcIlA\ngloEGJTgKg+4nQCDe4goKEyAAZCSAINynDCDGwkwuEosQRUCDB7n9Bg8QYDBbQQVlCTAAEhJgEEh\nTpXBvQQY3GAjnMwrwkMEGNxJXEEx+QKs/2PjDSXLA0AV/6ldgGP6vh+GYfn1+EqNQsFvToDB7fKN\nwDYMwzCLNChmjKjpLKIZRXhOUwG2qv9brTJU+dyZCMWIUIZOMYKVoVOMv1Vsqeo2lUeFnkKcbcRz\noytjMqozf0gK09YyRYaFDjDZQy5DN3zmDPuuH78GHhI6wJaGYRj7BRtXc0B10gue1njTn2IUTGtm\nh5R9kJzip0PjAQZVjMMvZ7/gOcmmECEFuQUFtH8ZPQBNEmAApCTAAEip5XNgywvug3x0mev+N8pQ\ncst8+6zCtZOlRla/W6YYEWrklrULLpZh9t33HCMZNRtg28v+1vroYpf1//zzP688vWW2i1GmDD+L\nUaZSdtZIrWJMG+uKNVLyYPm5V9Q9Riq2YImYQiwqyHLDyjCK0zREWIDuU4YIGyRIMQiu2REYP1Vv\nI6q316EUG/1sF6BuGUIpNoW4pxisEmBvFOTIHJvsWgX4fPT4/4obpHpdhFI9QSNM343L5k3Xz2PG\nFOJLVW8gKn76aPijq7pBgmwNQvkEp57NtpbnCuJchbi8lKPKpU2fz61+odesGFWueQtVI6vfKl+M\nIDVScsSjRrJrOcAAaJgpRABSEmAApCTAAEhJgAGQkgADICUBBkBKAgyAlAQYACkJMABSEmAApCTA\nAEhJgAGQkgADICUBBkBKAgyAlAQYACkJMABSEmAApCTAAEhJgAGQkgADICUBBkBKAgyAlAQYACkJ\nMABSEmAApCTAAEhJgAGQkgADICUBBkBKAgyAlAQYACkJMABSEmAApPT/ARADy1I9OKpHAAAAAElF\nTkSuQmCC\n", "output_type": "display_data"}], "prompt_number": 15, "cell_type": "code", "language": "python", "metadata": {}, "input": ["lw = 2; msB = 30;\n", "mystem = @(x,y, col)stem(x, y, [col '.--'], 'MarkerSize', msB, 'LineWidth', lw);\n", "subplot(3,1,1); hold on;\n", "mystem(x0(I), a0(I), 'k'); box on;\n", "plot(t1, y, 'LineWidth', lw); axis tight; title('Observations');\n", "subplot(3,1,2);\n", "plot(t1, [y-y0 y-y1], 'LineWidth', lw); axis tight; title('0th order');\n", "legend('0th order', '1st order');\n", "subplot(3,1,3);\n", "plot(t1, y-y1, 'g', 'LineWidth', lw); axis tight; title('1st order');"]}, {"source": ["Proximal Operators for C-BP \n", "----------------------------\n", "For simplicity, we re-write the initial C-BP optimization problem as\n", "follow.\n", "The C-BP constraint set reads\n", "$$  \\Cc =\\enscond{(a,b) \\in \\RR^N \\times \\RR^N}{ a \\geq 0 \\qandq \\abs{b} \\leq a}.\n", "$$\n", "We introduce the convex gauge $J$ such that\n", "$$ \\forall u = (a,b) \\in \\RR^{N} \\times \\RR^N, \\quad J(u) = \\norm{a}_1 + \\iota_{\\Cc}(u),\n", "$$\n", "and the linear operator $\\Gamma_z$ such that\n", "$$ \\forall u = (a,b) \\in \\RR^{N} \\times \\RR^N, \\quad \\Ga_z u = \\Phi_z a + \\frac{\\De}{2} \\Phi_z' b.\n", "$$\n", "Then, we may conveniently rewrite the C-BP problem as\n", "$$ \\umin{u \\in \\RR^{N} \\times \\RR^N} \\frac{1}{2}\\norm{y-\\Ga_z u}^2 + \\la J(u).\n", "$$\n", "\n", "\n", "The main ingredient to develop these schemes is the so-called proximal\n", "operator associated to a proper lower semi-continuous convex function\n", "$f : \\Hh \\rightarrow \\RR \\cup \\{+\\infty\\}$ (where $\\Hh$ is some Hilbert space)\n", "$$ \\forall u \\in \\Hh, \\quad \\text{Prox}_{\\ga f}(u) = \\uargmin{u' \\in \\Hh}  \\frac{1}{2}\\norm{u-u'}^2 + \\ga f(u').\n", "$$\n", "\n", "\n", "One can note that the $J$ functional, defined on\n", "$\\Hh=\\RR^N \\times \\RR^N$ can written as\n", "$$ J(a,b) = \\iota_{\\Cc}(a,b) + \\dotp{(a,b)}{\\xi}$\n", "$$\n", "$$ \\qwhereq \\xi = (1,\\ldots,1, 0,\\ldots,0) \\in \\RR^N \\times \\RR^N $$\n", "\n", "\n", "\n", "It is thus the perturbation of $\\iota_{\\Cc}$ with a linear form.\n", "This allows one to compute the proximal operator of the  $J$ functional as\n", "$$ \\forall u \\in \\Hh, \\quad \\text{Prox}_{\\ga J}(u) = \\text{Proj}_\\Cc\\pa{ u - \\la \\xi }  $$\n", "where $\\text{Proj}_\\Cc$ is the\n", "orthogonal projector on the cone $\\Cc$.\n", "\n", "\n", "It can be computed as\n", "$$ \\forall u = (a,b) \\in \\Hh, \\quad \\text{Proj}_\\Cc(u) = (\\tilde a,\\tilde b)\n", "$$\n", "$$ \\qwhereq  \\forall i=1,\\ldots,N, \\quad (\\tilde a_i,\\tilde b_i) = R^* ( R (a_i,b_i) )_+\n", "$$\n", "where $R \\in \\RR^{2 \\times 2}$ is the rotation of angle $\\pi/4$, and where we have noted\n", "$$ \\forall (\\al,\\be) \\in \\RR^2, \\quad (\\al,\\be)_+ = (\\max(\\al,0),\\max(\\be,0))\n", "$$"], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 16, "cell_type": "code", "language": "python", "metadata": {}, "input": ["C  = @(u)u(:,2)+1i*u(:,1);\n", "Ci = @(v)[imag(v), real(v)];\n", "ProjOct = @(v)max(real(v),0) + 1i*max(imag(v),0);\n", "ProjC = @(u)Ci(exp(1i*pi/4)*ProjOct(exp(-1i*pi/4)*C(u)));"]}, {"source": ["Callback for the proximal operator of $J$."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 17, "cell_type": "code", "language": "python", "metadata": {}, "input": ["ProxJ = @(w,lambda)ProjC( w-[lambda*ones(size(w,1),1) zeros(size(w,1),1)] );"]}, {"source": ["Numerical Resolution with Foward-Backward\n", "-----------------------------------------\n", "Starting with an initial $u^{(0)} \\in \\Hh$, the forward-backward iterations\n", "to solve C-BP reads\n", "$$ u^{(\\ell+1)} = \\text{Prox}_{\\tau \\la J}( u^{(\\ell)} - \\tau \\Ga_z^* ( \\Ga_z u^{(\\ell)} - y ) ).\n", "$$\n", "For a step size satisfying $0 < \\tau < 2/\\norm{\\Ga_z}^2$, one can show that\n", "$u^{(\\ell)} \\rightarrow u^\\star$ which is a solution of C-BP,\n", "see for instance [CombettesPesquet10](#biblio)."], "metadata": {}, "cell_type": "markdown"}, {"source": ["__Exercise 1__\n", "\n", "Compute the Lipschitz constant $L=\\norm{\\Ga_z}^2$. _Hint:_ you can use power iteration to estimate the largest eigenvalues\n", "of $\\Ga_z^*\\Ga_z$. We display here the evolution of the estimate of this eigenvalue with the number of power iterations."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAALkElEQVR4nO3d4VLq\nRgCAUdK57//K6Q8qzYWoiEDyZc+ZjqNU7WohH7tmwzTP8wkAav7ZegAA8AgBAyBJwABIEjAAkgQM\ngCQBAyBJwABIEjAAkgQMgCQBAyBJwABIEjAAkgQMgCQBAyBJwABIEjAAkgQMgCQBAyBJwABIEjAA\nkgQMgCQBAyBJwABIEjAAkgQMgCQBAyBJwABIEjAAkgQMgCQBAyDpz9YDABjRNG09gu/M89Yj+M7B\nAzbt/z4CA9n9EZGUgwfsdDrNr3wWMU3TS7//qxn/tp47fs/WGM3xAwY7sRaYWXXeYPVJgmc/337/\n/c+YBQw+NWZd7j8q1htAnYAxlkqTdAG+dfAnUJ4hHt4eguQuxvEkDp5mYAS8oVK7f6gC1wSMvXh6\npTQJjk3AeKunVEqZgJOA8SK/D5VKAV8TMJ7jsWKpFPAwAeNB9xdLpYBXEDDupVjArggYn1IsYM8E\njP8pFhAiYENTLKBLwAb1dbrkCtg/ARvLZ91SLCBHwIagW8DxCNjBraZLt4ADELBj0i3g8DIBmz4O\nyZeXqJkWB+l5nm8/YUCWCoFxZAJ2+ijT8mXWljFbfX8cplzAaDIBW23SedY1YK4udAsYViZgp7/X\nDM8uc7L7v+oYtbNUCDzd18fSHSoF7CpXd6boGMW6MOUCXmR5tEzE7J+tB3CX219l4pf7RNP03z9X\n5lm9gEE1ZmC3Jxne3nLIsxAtFQJ8phGw01qWrm45UrfOVudbAJxlAjaaq3pJF8AVAdujZb2kC2BV\n4ySOoagXwD3MwHbEsiHA/QRsL0y8AH7EEuIuqBfATwnY9tQL4AECtjH1AniMgG1JvQAeJmCbUS+A\n33AW4gacLg/wewL2biZeAE9hCfGt1AvgWQTsfdQL4IkE7E3UC+C5BOwd1Avg6QTs5dQL4BUE7LXU\nC+BFnEb/KtIF8FJmYC+hXgCvJmDPp14AbyBgT6ZeAO8hYM+kXgBvI2BPo14A7yRgz6FeAG/mNPrf\n8tooAJsQsF8x8QLYiiXE51AvgDcTsMddpl/qBfB+AgZAkoABkCRgACQJGABJAvZbzuAA2ISAAZAk\nYAAkCRgASQL2oKtLIALwZgIGQJKAAZAkYAAkCRgASQL2K3YxA2xFwABIEjAAkgQMgCQBe4RdzACb\nEzAAkgQMgCQBAyBJwABIErDH2cUMsCEBAyBJwABIEjAAkgQMgKQ/Ww/gXtPH1S/mj3Mnrm65/YSX\njeSl3x6Au2QCdlqEap7n89vz7ed0LT98dcMA2FxmCVGTAFgqzcCmhxbvrr5KCAFWPXaM3VApYFd/\n6/rRVwHwteXRMhGzxhLiDn+VsgiwrcYM7HzWxuX91VvedhYiAHvQCNhpLUtXt+gWwFAaS4gAcEXA\nfmZ/f4wDGJSAAZAkYAAkCRgASQIGQJKAPcIZ+wCbEzAAkgQMgCQBAyBJwH7ALmaA/RAwAJIEDIAk\nAQMgScAASBKwH7OLGWAPBAyAJAEDIEnAAEgSsHvZxQywKwIGQJKAAZAkYAAkCRgASQIGQJKA/YzL\ncADshIABkCRgACQJGABJAnYXl+EA2BsBAyBJwABIEjAAkgQMgCQB+wG7mAH2Q8AASBIwAJIEDIAk\nAfueXcwAOyRgACQJGABJAgZAkoABkCRg97KLGWBXBAyAJAEDIEnAAEgSsG/YxQywTwIGQJKAAZAk\nYAAkCRgASQIGQJKA3cVlOAD2RsAASBIwAJIEDICkP1sP4F7TxyUx5o+/R02Li2TM83z7Cc/4jz7r\nOwHwZJmAnT7KNE3TJVHLmK2+D8BRZQK22qTzrOvrXE1/T6O0DWDVVFt0ygTs7Gp2dZmTffEligVw\nj+XRMhGzTMBuJ1vKBDCy0lmIWz07EEqAHWrMwM65Wp5neHva4SvOQgRgtxoBW23S1Y26BTCU0hIi\nAFwI2KcK5+AAjEvAAEgSMACSBAyAJAEDIEnAvuHkfIB9EjAAkgQMgCQBAyBJwNbZxQywcwIGQJKA\nAZAkYAAkCRgASQL2FbuYAXZLwABIEjAAkgQMgCQBW2EXM8D+CRgASQIGQJKAAZAkYAAkCRgASQL2\nKZfhANgzAQMgScAASBIwAJIE7JrLcAAkCBgASQIGQJKAAZAkYAAkCdg6u5gBdk7AAEgSMACSBAyA\nJAH7i13MABUCBkCSgAGQJGAAJAkYAEkCtsIuZoD9EzAAkgQMgCQBAyBJwP5nFzNAiIABkCRgACQJ\nGABJAgZAkoABkCRg11yGAyBBwABIEjAAkv5sPYB7TR/bjOePNb6rW24/AYADywTstAjVPM/nt+fb\nz+lafvhAw1yGA6Als4RoXgXAUmkGdnpodjX9PbcSQoBVU20lKhOwq3XC+ykWwD2WR8tEzDJLiCcp\nAmChMQM7PxdYnmd4Po/j8uHpSWchSiRARSNgq026utH8DGAopSVEALgQMACSBOx0sosZIEjAAEgS\nMACSBAyAJAEDIEnA/mcjGUCIgAGQJGAAJAkYAEkCZhczQJKAAZAkYAAkCRgASQIGQJKA/ccuZoAW\nAQMgScAASBIwAJIEDICk0QPmMhwAUaMHDIAoAQMgScAASBIwAJIE7HRyGQ6AIAEDIEnAAEgSMACS\nhg6YXcwAXUMHDIAuAQMgScAASBIwAJIEzC5mgCQBAyBJwABIEjAAksYNmF3MAGnjBgyANAEDIEnA\nAEgSMACSRg+YXcwAUaMHDIAoAQMgScAASBo0YHYxA9QNGjAA6gQMgCQBAyBJwABIEjAAkoYOmMtw\nAHQNHTAAuv5sPYCfmaZp/pg3TYvNXPM8Xz6cTawABpAJ2LS293gZs9X3ATiqzBLiPM+3WZqmaTVs\nX3MZDoADyMzAVp2T9nXDrv6tyRnAqgfmA9sKB+zOFCkWwD2WR8tEzDJLiFcSv1wAXqc6A7s97dBZ\niABDiQVsGaerUP20WzIHkFZdQgRgcAIGQJKAAZA0XMCcvQhwDMMFDIBjEDAAkgQMgCQBAyBp0IDZ\nxQxQN2jAAKgTMACSBAyApLECZhczwGGMFTAADkPAAEgSMACSBAyAJAEDIGnEgLkMB8ABjBgwAA5A\nwABIEjAAkgYKmMtwABzJQAED4EgEDIAkAQMgScAASBouYHYxAxzDcAED4BgEDIAkAQMgaZSA2cUM\ncDCjBAyAgxEwAJIEDIAkAQMgaayA2cUMcBhjBQyAwxAwAJIEDICkIQJmFzPA8QwRMACOR8AASBIw\nAJIEDICkgQJmFzPAkQwUMACORMAASBIwAJIEDICk4wfMZTgADun4AQPgkAQMgCQBAyBJwABIGiVg\nLsMBcDB/th7Az0zTNH+0aPo4v/B8y9WHABxbJmDT36fD35Zs+aGGARxeZglxnmdZAuAiMwN71Pp6\nIwBXptp1Hw4fsP/oFsDXlsfJRMwyS4gAsFQN2DzP04fzn8eWH75tGIknKV8w/m0Z/4bSgz/1x/8U\nsSXEZZyuQmWREGAo1RnYj0gbwPEMETAAjufge34/VoktFgP8zP7rcPCAAXBUlhABSBIwAJIEDIAk\nAQMgScAed7n2x9YD+ZXu+Ou///rgL+/kfpDlaOvj/+yWQcSuxLE3l9fSjJ7M2b3fL18Ervj7v4w5\nN/iro3/rdfhW7/Ch/xGr4+8+in/PDOxx+7+7fy3xiP3a+69+SfqV+W4H3/pZbsc/+P3fDOy3Br8D\nbaj13P/W1SuJs6HoXQgBe1z6AHQe/OVt9KeIyq28HZhHcZqA/Ur3HuMACmfd+79HsYA9aPnc51R+\nDESdXwHu8v62g/mp9OCX6j+IR3HdoN0GoM5ZiAAkCRgASQIGQJKAAZAkYAAkCRgASQIGQJKAAZAk\nYAAkCRgASQIGQJKAAZAkYAAkCRgASQIGQJKAAZAkYAAkCRgASQIGQJKAAZAkYAAkCRgASQIGQJKA\nAZAkYAAkCRgASQIGQJKAAZAkYAAkCRgASQIGQJKAAZAkYAAkCRgASf8Cvb343E+N52AAAAAASUVO\nRK5CYII=\n", "output_type": "display_data"}], "prompt_number": 18, "cell_type": "code", "language": "python", "metadata": {}, "input": ["exo1()"]}, {"collapsed": false, "outputs": [], "prompt_number": 19, "cell_type": "code", "language": "python", "metadata": {}, "input": ["%% Insert your code here."]}, {"source": ["Regularization parameter $\\la$."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 20, "cell_type": "code", "language": "python", "metadata": {}, "input": ["lambda = 40;"]}, {"source": ["Function $F(u) = \\frac{1}{2}\\norm{\\Ga_z u - y}^2$\n", "and its gradient $ \\nabla F(u) = \\Ga_z^*( \\Ga_z u - y ). $"], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 21, "cell_type": "code", "language": "python", "metadata": {}, "input": ["F = @(u)1/2*norm(y-Gamma(u), 'fro')^2;\n", "gradF = @(u)GammaS(Gamma(u)-y);"]}, {"source": ["Energy without the constraint $ E(u) = F(u) + \\la \\norm{u}_1 $."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 22, "cell_type": "code", "language": "python", "metadata": {}, "input": ["E = @(u)F(u) + lambda*norm(u(:,1),1);"]}, {"source": ["Initialization $u^{(0)}$ of the algorithm."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 23, "cell_type": "code", "language": "python", "metadata": {}, "input": ["u = zeros(N,2);"]}, {"source": ["__Exercise 2__\n", "\n", "Implement the forward-backward algorithm. Monitor the decay of the energy\n", "$ \\log_{10}( E(u^{(\\ell)})/ E^\\star  -1   ) $\n", "where $E^\\star$ is an estimate of the minimum energy (obtained by\n", "running the algorithm with a lot of iterations."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAANA0lEQVR4nO3d7XKb\nOACGUdjJ/d+y9odbSsDfxkivdM50utl0Z60S4icSsplLKRMApPmv9gAA4B0CBkAkAQMgkoABEEnA\nAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCS\ngAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIv3UHsB3zfOff9Yd\nBkCcUkrtITzQecAu6n4Z5nlu4TwwjNaG0cIYDKPBYbQwhsswag/hMUuIAEQSMAAiCRgAkZpYbP2e\nyypu139FgOM1cinuPjMwACIJGACRBAyASAIGQCQBAyCSgAEQaYiAJbwlCgCvGSJgAPRHwACIJGAA\nRBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgUvcBcy9L\ngD51HzAA+iRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEGmUgM1e0AzQl1ECBkBnfmoP4Ejz33lW\nKaXuSAD4tt5mYKWUUspsxRCgd10FzMQLYBxdLSFOq1XEW38kcgBXxa1dpQZsc6CXLF0+uPplkC6A\nO9ZPkhExSw3YvkbzPEsUwDhSA7a33ruhZADd6ydg041uleJVzAAd6moXIgDjEDAAIgkYAJEEDIBI\nAgZAJAEDIJKAARBpoIB5NRhATwYKGAA9ETAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACINETA3J8Z\noD9DBAyA/ggYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBorYG7KDNCNsQIGQDcEDIBIAgZAJAED\nIJKAARBJwACIJGAARBIwACIJGACRRglYKbVHAMChRgkYAJ0RMAAiCRgAkQQMgEgCBkAkAQMgUocB\nm922EmAAvQXsYb3UDaAPXQVsnufiFcsAY/ipPYCvW83JtA3gprjrL6kB2xzoUsrlM8vvy1Rs+SDt\nSwNwqvUKVkTMUgO2XypchcpCIkD/uroGBsA4OgyY6RfACDoMGAAjEDAAIgkYAJEGCphLYwA9GShg\nAPREwACIJGAARBIwACIJGACRBAyASAIGQCQBAyDSiAFLuM0NAA+MGDAAOiBgAEQSMAAiCRgAkQQM\ngEgCBkAkAQMgkoABEGmsgLkpM0A3xgoYAN0QMAAiCRgAkQQMgEgCBkAkAQMgkoABEGnQgLmnJUC6\nQQMGQDoBAyCSgAEQScAAiCRgAEQSMAAiDRcwd1QB6MNwAQOgDwIGQCQBAyCSgAEQadyAeTtEgGjj\nBgyAaD+1B3CkeTWrKvbLA3Stq4BNugUwjN4CdpmErTM2/77YpXAAV81pWwN6C9ilT/M8L6FSLIBn\n3PnRv02pAbs6r3qyVaXYgggQLzVg+1atZ10AdC81YHullGVapmQA3esnYJNuAYzEC5kBiDR0wGzl\nAMg1dMAAyCVgAEQSMAAiCRgAkQYNmP32AOkGDRgA6UYPmJ30AKFGDxgAoQQMgEjjBmzZx2EVESDR\nuAEDINrQATMJA8g1dMAAyDV6wEzCAEKNHrBJwwAyCRgAkQRsmkzCAAIJ2B8aBpBFwK7QMID2Cdg/\n63usaBhA4wTsFw0DSCFgWxoGEEHArtAwgPYJ2HWbhskYQGsE7KZ1wyZTMYDGCNg9pVhOBGiUgD1m\nORGgQQL2FMuJAK0RsGdZTgRoioC9xnIiQCME7GWWEwFaIGDvsJwIUJ2AvU/DACoSsI9oGEAtAvYp\n2zoAqhCwA7gkBnA+ATuMhgGcScCOpGEApxGwg2kYwDkE7HgaBnCCrgI2/1V7IBoG8HX9BOzSrVJK\nKUXDALr3U3sAB1syVnsg0/S3YZd6zfP2TRQB+ERXAVu6Nc/z+uOr/82Jo/rXsGn3XsAAjWhh7eol\nqQF7PkstzMaWhk2mYkCr1s+WETFLDVgLWXrJpmGTqRjAZ1IDtrfeu9Fm3taXxCYZA/hMPwGbWu3W\nxnoqNskYwLv62UYfZPPmv5O3sQd4nYBVI2MAnxCwymQM4D0C1gQZA3iVgDVExgCeJ2DNkTGAZwhY\no2QM4D4Ba5qMAdwiYAFkDGBPwGLIGMCagIWRMYALAYskYwACFuxqxgAGIWDxNhkzFQMGIWCdsKII\njEbA+uHCGDAUAeuNC2PAIASsTy6MAd0TsJ7JGNAxAeufC2NAlwRsCC6MAf0RsIFYUQR6ImDDkTGg\nDwI2KCuKQDoBG5epGBBNwEZnKgaEEjBMxYBIAsYfMgZkETB+saIIpBAwtkzFgAgCxnWmYkDjBIyb\n9lMxgHYIGA9YTgTaJGA8ZioGNEjAeJapGNAUAeMFdnYA7RAwXmM5EWiEgPEOy4lAdQLGm0zFgLoE\njI9oGFCLgPEpy4lAFQLGAexOBM73U3sAh5l/P2uWzXMqX3Y53ssXYZ63VQM4Vj8BWxdrNgWopBQN\nA07S4RLiPM+mXxXZ1gGco58Z2C2WFs+3mYdNu4tkQIPi1q5SA3YrS/vpl2JV4ZIYxIm7EJMaMFmK\n4JIY8D1dXQNz9atBLokBX9JVwNSrTRoGfENXAaNZGgYcTsA4iYYBxxIwzqNhwIEEjFOtb8KiYcAn\nBIwKNAz4nIBRh4YBHxIwqtEw4BMCRk0aBrxNwKhMw4D3CBj1aRjwBgGjCRoGvErAaIWGAS8RMFqk\nYcBDAkZDvNcU8DwBoy0aBjxJwGiOhgHPEDBapGHAQwJGozQMuE/AaNe6YQAbAkbTvDgMuEXAiKFh\nwJqA0ToXw4CrBIwALoYBewJGBhfDgA0BI4aGAWsCRiQNAwSMJC6GAQsBI4yFROBCwMijYcAkYKTT\nMBiWgBHJxTBAwEhlIREGJ2AARBIwgpmEwcgEjGwaBsMSMAAiCRjxTMJgTAJGVzQMxiFg9MDLwmBA\nAkYnLCTCaAQMgEgCRj9MwmAovQVsnufZsxfAALoK2DzPpZRSioYNyyQMxtFVwGBNw6BvP7UHcLDL\n3KusdlVvZmPFhuvelfIvXfNshz08K27tKjVg+yxd1g+XP10+VqwBbRq2/jxwy50f/duUGjBZ4r51\nwxabzziJIFpqwPbWezfkjemJDR16BtH6CdikW9ywOS/0DPrQVcDgGXoGfRAwRqdnEErA4Jd1me7s\nw7K5EaoTMLjJ5AxaJmDwLD2DpggYvEnPoC4Bg2PoGZxMwOAr9Ay+TcDgDHoGhxMwqOCNnokZbAgY\n1PfMi89MzmBDwKAtFhvhSQIGTdMzuEXAIImewULAIJjNIIxMwKAfz/TM5IxuCBh0S8/om4DBKGzW\npzMCBiOyGYQOCBigZ0QSMGBLz4ggYMADNuvTJgEDXmNzI40QMOAjekYtAgYcyWZ9TiNgwLfYDMJX\nCRhwEj3jWAIG1KFnfEjAgCbYrM+rBAxokc2NPCRgQACbG9kTMCCMi2dcCBiQTc+GJWBAV/RsHAIG\n9EzPOiZgwEBs1u+JgAHjsrkxmoABTJPFxkACBnCFnrVPwAAe07MGCRjAy/SsBQIG8CmbG6sQMICD\n2dx4jq4CNv89I4pzAWiDxcbv6Sdg8zwv3Vp/DNAOPTtQPwG7Zf59Iggb0I6mejbfevhWdRgwxQJC\n1e3Z+tkyImapAdtXqpRy+eTyAUA0mxvvSw3Yfl7luhfQN5sbN/6rPYDDXCZeF02VrJHpoGGstTCM\nFsYwGcZvLQzjyTGU8uvX7f/br1+dSZ2BXdVUtwBO09RmkNN0FTAApmF6JmAAnXurZwE1a+ty0eFa\nWNEGaNv1CrQfh85nYH3nGeBYWT/zdx4wAJ6X9TN/P9voARiKgAEQScAAiOQa2Fesdz9eNpKcf6+y\nzf1lag1mGUatY/Lw737CMO486J2BnTmMFr4iJw9jGUzFE2Mzhqn2iXHrQZu91aKAfcvmfZ3PvFfZ\n5sUD+2+Pcwazfw1DrWOy/ibcPOhpw1jGsDxElaOxHkatE2PaHY1aJ8Zyiu4f9LRhNPJtcudBpxNP\njFdZQvyWy7syVnnoy3vzV3no+8OockwaORT7T55/NFo4FNONd+I+/8Ro4bn46hhqHY3EV82agX3L\n/iduKh6T5cfqkx93P4blX6scjUaepK7etM83y0WVoxH6JRCwr8g6Cc5R65hs1kAaGUOt8ewvO1Uf\nRpVDcXno9e/Vx7BfTT1N7vOVJcTjVX92aFDdY9LC92f1e902clruL9BWGUb5a6r6w8RmDCOfGO8J\nmzCmqLijaT2GNnchVtn+tzzQ+cO4tVx28jCeedBx9oWuB1N9GHe+Vc8ZRiMnxhsEDIBIlhABiCRg\nAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJ\nwACIJGAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMg\nkoABEEnAAIgkYABE+h+DjUE5isAH6AAAAABJRU5ErkJggg==\n", "output_type": "display_data"}], "prompt_number": 24, "cell_type": "code", "language": "python", "metadata": {}, "input": ["exo2()"]}, {"collapsed": false, "outputs": [], "prompt_number": 25, "cell_type": "code", "language": "python", "metadata": {}, "input": ["%% Insert your code here."]}, {"source": ["Recovered parameter $a,b,x$ that are intended to approximate\n", "$a_0,b_0,x_0$."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [], "prompt_number": 26, "cell_type": "code", "language": "python", "metadata": {}, "input": ["a = u(:,1); \n", "b = u(:,2);\n", "delta = Delta/2 * b./a; delta(a<1e-9) = 0;\n", "x = (0:N-1)'/N + delta;"]}, {"source": ["Display recovered measures."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAATwUlEQVR4nO3d25bj\nKLYFUOuM/v9f1nlQpsrpWzhsIfaCOR+6o+zqDsIIFiCEl3VdLwCQ5v96FwAAPiHAAIgkwACIJMAA\niCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACI\nJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIj0v94FgMEty7L9sK5r35LAYBaNChrZo+ua\nFgdHGXAJ8WGvASd7dh26PuEoQy0h6hoA5jHUDGxdV+szVPB6LGWkBYcYagZ2T09BTa5M6qs/Hxg8\nwC416mBZSmyWqVCMCmU4oRg/5tP22yt8GhXKoBjVynAJGWMNtYQIRbzugCp0TzAAAQZApAEDzPCW\nCp5dh65POEqJxdZ2iqwmM6e/NxH2ewnr5XJxPRIhovMcfxMH9LX3Agk3xSHJgEuIUIG4gtYEGJxK\nsMFRBBgAkQQYAJEEGDS0beBYlmVZlvJbuiCMAAMgkgADIJIAg+O93mpoIyIcQoABEEmAARDJUVLQ\nXP0z5SCRGRi0ch9bggwOJMAAiCTAAIgkwACIJMDgYO885uVRMPieAIPmtrMQe5cCRiPAAIgkwACI\nJMAAiCTAoIlnzyx7lhmOIsAAiOQsRGjOWYjQghkYAJEEGBzp/ce9PBgGX8pbQtwfCL1flnnxFgCD\nCQuwZVn2cLr++fVbAIzHEiIAkcJmYO+4OXTu5h/NzDjfdhG69igu7sTOcQJsXde9m7iuBr0GwDtu\nbsp0LMmbxgkw972o4/WVuK62IMIBwgLsenZ1s2Xj4VsAjCoswC6Pwml/RW4BzMMuRAAi5c3AII61\nAWjBDAwO89utGbZywDcEGACRBBgAkQQYAJEEGACRBBg0tyxLxME8kEWAwcHsmYdzCDDoQMjB9wQY\nAJEEGACRBBgAkZyFCM05CxFaMAODY3y2T97ueviYAAMgkgADIJIAAyCSAAMgkgCD5pyFCC0IMAAi\nCTA4kie+4DQCDPoQdfAlAQZAJAEGQCRnIUJzzkKEFgQYwOX+OQfDjvosIQKze/iUnkf36ssLsOWv\nF++eXCT45qJzwZalMykubAlxWZZ9Xn/9880r928BPCSlcoUF2Adurk7BBrxvqtFwXJaPFmBbBVxf\ncPNcfJR1f1mSYqpau/5jI8JsnAB7vboIwGDyNnEAHMhgN1fYDGxd131iezPfevgWwMf0JMWFBdjl\n0SW1v+JqAz6wruuyXC6Xm7s++pPqLCHCYX47gjLiKmJfu/nzX+sqvSLkzcAgjrWBFHePlhpklGYG\nBkAkAQZwuVjRDSTAgKk9fGBXmEUQYABEEmDwre/P3Ek4tQfKEWDQnG/5SXFfU+qtMgEG4KZXJAEG\nQCQBBszrxQqhOVl9AgyASAIMgEjOQoRjvFhxchZiioc15UTEsszAAIgkwKAnQ/sK1EIoAQZAJAEG\n8JiZWXECDJiUY6LSCTBozlmIKdRUFgEGQCQBBl85arxu3F+Z2qlJgAEQSYABU7PVMJcAAyCSsxCh\nOWchplBTWczAAJ6SaJXlzcD2pzRuxko3T28YSQEv2Fg4gLAAW5ZlT6brny//JpZHEQGGN+AS4k2w\nwQlccXC+sBnYBywtUty6Ws4KMMPXWsatXY0WYPfTL4lFd1u/4FKsb/KairsRM+ASIsCbZo2qQYTN\nwNZ1vd+FuM+63P0CmEdYgF0eze73V6QXwDwsIQIQSYDB5xLuc/PAryrOyk5ZeUuIEOfNxe0ZNmoX\n5zZEFjMwACIJMIC3WDGuRoABEEmAARBJgAGTsmMjnQCD5pZliThZDjWVRYABEEmAwbe+X4mylgUf\nEGAAPzDCqEmAAXNxk2sYAgyASM5ChOacsJdCTWUxAwN4l+XHUgQYAJEEGHzIYBz6EmBQiFA8jbtd\nAxBgAEQSYNCcE/ZSqKksAgyASAIM4GfumRUkwICJWCAciQADIJKjpHq6uV3sGBuA9wmwbu43Oy3L\nIsPivFNj71Trulrd6k8DzJIXYHu/f3+pvXirmmdbdWUYFLcsNnRUERZg1/37TV+/RcL2SvEY8KAJ\nwPdG28SxPYdYOb1+JN6gteQegv+EzcBeezg5s1EC4B1xo+ehAuyhuMSKKzAwhvubMsWNtoQYQUTN\nxgl7KdRUlrAZ2Lqu91sNtwXDh28BMKqwALs8Cqf9lajcWi+XhwO9oD9hau2G6XZpl+VZvWosIXa0\n7nG1ruv1PwItiJ/BCDAAIgmwnq5XiqwaAfxK3j2wAVyvYzw6EEuYjSbq7uzU1FQWMzCA33EvrQgB\nBkAkAdbN/VqF1YvJuQDgVwQYfE7kJFJrw7CJ42zvrJ7bxzGUfzft9CsHjMYMrDNnrw3upnKXxQaA\nyrTHLAIMmnnWFeoi4QgCrA8rSeOTUiPScksRYLVoHgBvEmCnMiiHXrS+8QiwojS2ytQOVGAbfWfO\nXuOe5yh60R6zmIFBGy+6Qr0kHEGAdaD7moWaHpdl5AosIZbja8uHsmXYVqPyDA4lwKA90VWGqhiJ\nJcTzmFcBHEiAdfbi7DWBNw5HIIZwFmIWAQYfarEYZYEL3ifAAH7BIKMOAXY2Vz/AIQRYRUIOjuXG\n1pAEGACR8p4D2/cI3Z9adr19qNqZZs8GgNXKSRNqOYT2mCUswJZl2a+w6593g11/DnUFeGa0JUSP\ncQDn0NN0FzYD+9E2A7uZqN3/C/Cxc7otk2/OFzf6HyrAHoaTxAJ4x3VvGRFm4ywhRnzcwhQ60gAH\nEzYDW9f1fhfitmD48K36tjI/2o1ihX0gvk4lxLP2SE1hAXZ53Nevz94CYFTjLCFWZi4FcDgBVp3w\ng2qs9RQhwACIJMDgE+3G4Eb3h7OMMaq8TRyDsfFkCmo5hPaYxQzsPJoGwIEEWF0CD+AFAQbwIXfX\n+hJgAEQSYM0ZowG0IMA6e+cLzERgHR/WxbJ88L9U78d656ayLxTMIsAAiCTAAIgkwACIJMBO4qEu\nGIkWXYEAK00jgS/ZkzEwZyF25uy1KajlENpjFjOwtoz+ABoRYPBrrYfppgFn2x7UM95MI8AyaFnQ\nxE1uibEoAgyY1bOs+k2GybuOBBjAHbmUQIB15uy1KViY6s1txSEJsDNoPFCOIUU+AVad8IMmNK18\nAgzedf6Q3SQBXhBgAHfMzxIIsIYMn6G0Zyn1XnrJuO7yzkLc9+w9O7VsWZagA83eL+qyaDCx1Fwn\nPw8it6rZ/z01FSUswK7D6WFQ2ZIO/JrcyhQWYK9tkXaTYTf/GDQ5AzhT3ARgqAB7SGIBvOO6t4wI\ns3E2cWwf9/V/FiFAAVoYJ8DWvy7DzbrG+mvinVMdKj1IpQHzXMKWEK9vcb3ezZFi+3Nyy89btotW\nLVe2Ld5cLhftMUdYgF0eXVs3r7j4gGu6hFGNs4QIwFQEWCuWxQGaEmBJhCLAToDBW4weuOfuWl95\nmzgGY8vJFL6oZWdgfux3Y47tCZxGRaENMzAAIgmwtoydARoRYBkEIcANAQZAJAEG8C2bVLuwC7GJ\n969mZyFO4aOzENdVt3giZyEGMgMLo0frS88WR5UNTIABEEmAARBJgAEQSYABfM49to7sQmzonSvb\nfqcI3+6d+a6WHYf4gV9XmbMQA5mBxdCFAVwTYABEEmDH86gWwAkEWB4BCQVpmOcTYABEsguxM2ch\nTuGjsxD53i8+cmchBjIDg3fp1qAUAQZ1iUx4QYC1ousBaEqAJRGKALu8TRzL382q9zdaX7wF0Iiv\nHu0lLMCWZdnD6frn3fbKw7dqSiknX1HL5/okTpyFGGioJcQKYWAgNp4KdVqhDFBN2AzsR8tdQ795\npULIfc/x5MDh7vvP4kYLsH0J8eYVAF677i0jwmycJcSIjxsYmE7oZGEzsHVd77cabls2Hr4FwKjC\nAuzyKJz2V+rk1vsFcRbiFJyF2MPvPm9nIQYaZwlxEloWwEaAwVt6DR0MWeAZAQZAJAF2JHuQYE4m\nyl0IsFTCEh7SNOaRtwtxMPY7TUEt1+csxEBmYPBKneF8nZJAEQIM4DDGGWcSYMezXARwAgGWR0AC\nXAQYMCTjvBnYhdiZsxCn4CzE+pyFGMgM7DBu3gKcSYAFE5mn6TsiNx+IoJrOZwmxk/2ry/oWA8Zi\nVDcVM7Ae7huZZgfwSwLsdM+ySoYB/IYlxIP9sA4upaIcVl0H3R5ZFjdamjnuLETVdBozsEre7i81\nDwABBkAkAVaJiRUcQUuahAA7xruLfxoWwEEEWDabQqAOA9STCbDTPbvGXftVHVAzy/LlWMPV8Y6v\nPuNluSx/HFYgGhNgPazr3iEt2/mh+ieAXxJgRxJDAKcRYKmEJTC5vJM49hXq++/sefEW/FbNWyFO\neYigms4RFmDLsuzhdP3zbnvl4VvADDT9eYQF2Gu9QuvDofpBZ68Z6wVQQ/UddxYipxkqwDY306+b\nTbFmZgAPxT1CMFSAbZ/+TURJLIB3vBj61zTaLkRxxbFcUCmK9LcumDOFzcDWdb3faritGW6vd9yI\n6MKlqXWt0kdDEWEBdnmUTNsrE8699GjAzPICbBzblPFyuUyZvnPZBhpqubKj26PtwScY7R7Y+cyB\nhqRaoT4BNgK97VRU9wsmPVMRYABEEmAARBJgAEeyjHkauxCP8ckle8TZa3bSN3VYT6RLa+yAVuAs\nxEBmYBBDDsI1AQbQhNWR1gTYV+pcoHVKAnAOAQa36o8G6pewC0ussxFgAESyC/EAH477nIU4D2ch\n1qc9BjIDi6etMbmCC6pa5TkEGDymD4LiBNjnCo77GJ5YzaKXaEqAjUNTAaYiwIARmJtOyC7Efpy9\nVlKTiWyDztUX/h5MewxkBvatCp1IhTIA17TKEwgwIJhbvzMTYB/SbMZWefhcuWxwJgE2FLEK1WiV\n7QgwACLZhdiPs9fqaTVYbnMWoo2IuwM+B+0xkBnYV+pc53VKAnAOAQakKn57ybCytbwAW/568S+0\nL0Pr3/C5ymUDOFDYPbBlWfbl6euf91d6FIrR1B84r6uRShJ3KxvJm4G9sK6ru68AkwibgX3gZlp2\nYMJ9+//k7LViGs5pmo2rDO0vR3262mPgItb4ATbPnMyyEvCN694yIsyGWkI8R/1qrV9C+F7EdT7N\n+LmPsBnYuq77uODFbg6AUiz2thAWYJdHS4I3rwgzvpFy+VgxBkuIH6rZzdUsFbTjmp9Z3gysryPH\nvM5eq6TtbKbNWYj7/7fL5wDN2qO5cjtmYGPSYBibK5yLAIPz6HTnpv4PJ8A+YcVmVMfX7LL802/d\n/ON3XIc+gckJsF+IGEBp0h9o+zVg7UVcmXA4ATYsnVpnLypA3Xwn7vMzrGzELsR+nL0GdZzSHm0Z\nPZYZ2K/Vv/7ql7CmxM8tscxHmflvZyPA3hW3anHJLPP5xviUxvgr4FcEGJAkNKrNF1sQYNDGix5L\nZzax0ACuSYC9Zb/mUnqelHLW0eQTO6Ua5qzrOf9qbtiF2M8pZyHa9fRa8+Hw9umfMgKaoa4b1lf7\n9uhQxMOZgcFZho8X3iPGjiLAfha3frjZS6u1/CirZu+ll/+3Zvt7eUaAMa/xon28v+jaAH+d6D2W\nAPtB6PRrk1hmmMEAYVyBAJuC1nIvemhyb/gV48Hqi0PYhfiWdnusNcYp6HTrO6s97nsRZ9g12poZ\n2CsDDGaHH5h/Zsjh/MB1PWR98T0B9tR4vQBQwcBDjZMJsJ+lD/q0lhsDD+eHrOuB6+syVk2dT4A9\nNlibGbJfYwajXrFjdCzdCbAHRm0zm7H/uh8NNjS5N+pgZbz6GrWmziTAbl1fTG3bzLJclj9a/po/\nrv+WaRtMt/RaljM/9GF6xlPr69z2eP/L+YAA+0ej9OrSJB75rxi9StTxo7j+zUVq5IRi/Pgban4U\n14k/3txr93BYWaRGIgiwP26GyKO2mZsGM0lLmaRyr0VX9FSVFV1T3Q31IPM+cnn/2xAeXi5jt5mb\n73TYfh7yT76v3CH/zGduvsil/oRm2vq6+5qV1TPObxonwJZl2XPr+ufL3yfe3xnaTHLR3HRtlyd9\nR0QrUrmv3X8495/VmRX9vL5uSzBVfY3UJM+0tPsqxZM9DLDKq8lbWeuWj+Oo6/rU0b366TDODOyJ\npxVQpGpqlOKPwnH/a0Xq91qpEpWq6zqVVaYg/ylVU9UMHmB1GkYEH9c81HUKNfWCXYgARBrnHtjl\no12IAIQaKsAAmMfI98A6Tshe/+qbXf7nl+HMT+bZ7zq5dlJq5OG75xSjQo3cbBvu9VFcvztPG7l3\nTrv4xrAB9uKxsI6/+rRt/T/++ftjBk0/mdfFOKcMPxbjnEp5s0Z6FeO6s+5YI2c2lh+vir5tpGMP\ndqn9ANI1mzhOta5rhRGNMuzqjDF7HSN7X4YKH0iRYkyrSE/1o2FnYPyoex/Rvb8u5bTZz+sC9C1D\nKactIb5TDB4SYDMq0jK7n5ay/er9Pzt+IN3ropTuCdp3+W6zruveTmXYM5YQJ9W9g+j423frX5eu\nH0iRT4NStuA0snlt5LWCOrsQ77dydNnatP3evhu97o+p7LLnrVSNPHzrtGJUq5EzZzxq5MeCFQ+I\n6uUDgIcsIQIQSYABEEmAARBJgAEQSYABEEmAARBJgAEQSYABEEmAARBJgAEQSYABEEmAARBJgAEQ\nSYABEEmAARBJgAEQSYABEEmAARBJgAEQSYABEEmAARBJgAEQSYABEEmAARBJgAEQSYABEEmAARBJ\ngAEQSYABEEmAARBJgAEQSYABEEmAARBJgAEQ6f8BwnxmIPkzaYUAAAAASUVORK5CYII=\n", "output_type": "display_data"}], "prompt_number": 27, "cell_type": "code", "language": "python", "metadata": {}, "input": ["J = find(a>1e-3);\n", "t = (0:N-1)'/N;\n", "s = (0:P-1)'/P;\n", "clf; hold on;\n", "plot(s, y, 'LineWidth', 2);\n", "mystem(x0(I), a0(I), 'k'); % initial spikes \n", "mystem(t(J) + delta(J), a(J), 'r');  % recovered spikes\n", "axis([0 1 0 1]);\n", "box on;"]}, {"source": ["Display estimated vectors."], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAASzklEQVR4nO3d25Kj\nvBkFUEjl/V+ZXDDjMIBpcHPQ/rxWpVIeN61famE2EjL0wzB0AJDmP09XAAA+IcAAiCTAAIgkwACI\nJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgk\nwACI9N+nKwDfqO/71+thGB6sCeQSYPCAV2hNkww4RIDBA+QW/J4Ag7v1fW8EBr8nwOABcgt+r3cB\nGYBEltEDEMkUIlzLinm4iClEuMq7C10+dHCKglOILo8DfINSAdb3vfSiERu7or0UTlEqwIZhMDkD\n8CWKL+Jwqkub7Jm0r/3xQPEA6xL64GV6g4ZKvrNd2xHV8h/kO/sr10XtijjHKjWFCMD3EGBwiY2T\n4pLjALhfwQBzdKARq7ui/RPOUnNS+KXqrDdZxssJdkWCRBw86y/igMe1fyCARAWnEAH4BgIMgEgC\nDIBIAgwu5y6dcAUBBkAkAQZAJAEGQCQBBkAkAQZAJHfigMu5EwdcwQgMgEgCDIBIAgyASAIMgEgC\nDIBIAgwu516IcAUBBkAkAQZAJF9kBvjEbFrY19XvZwQGcMzqRU2XOe+XF2D9Xxs/vblKAJ0Mu13Y\nFGLf969x+vT17J3lj+BB9sZKpFQ7wgLsA+apgdtEnz3HZXO1ABs7YLoD5e5MAHeaHi0jwqxOgG3P\nLgLcwJHnTnmLOAAeJKLaETYCG4bhNbCdjbdWfwRwG0eem4UFWLe2i7zesffQpuWlWaK9Vjsv3+RO\nphABPjE9dZZejxBgAEQSYABEEmAARBJgAETKW4UIcVzhr0rPPssIDIBIAgyASAIMgEgCDIBIAgyA\nSAIMLtf3fcTTlThKzz5LgAEQSYABEEmAARBJgAEQSYABEMm9EOFy7phXlZ59lhEYAJGMwFrx+jaJ\nczqAPQTY82ZfhBz/KcYAtuUF2MZIJXEQ8+5r/H3fB7UC4H5hATY9rM8O8dOBi6M/QHnVFnGMtyZL\nSa/tu6i5x1oZ7phXlZ59VtgIbNvq4Gy2e6VkG8DN4sK4VICtklgAeywvyjSu2hRilu1wFb0AG8JG\nYMMwLJcajhOGqz8CoKqwAOvWwun1TmJuTXN39v79leEi+hKukBdg9byW/k//SRHTs5Pxtf4txKf1\nWQKsFT4JAIdYxAGXWV3HlbC4CyIIMAAiCTC4zOq0sLliOIkAAyCSRRwN8SCVgiZLTPVsPXr2WUZg\nAEQSYABEEmAARBJgAEQSYABEsgqxIdYyVaVnq9KzzzICAyCSAAMgkgADIJIAAyCSAAMgklWIDXFf\ntar0bFmvp7vp3CcYgQF8ZPpsUs8pfYIAAzhumVgy7HYCDIBIeQHW/7WxwZ31AeARYYs4+r5/XQmf\nvp5ucHulAHhAWIBtGyNtlmGzf7a8EqzluvEberagYZhf9Mrv5bgBQKkAW+XYAVyi3LFlerSMCLO8\na2DvjH/u6f8DUFidEdj2tTEAigkLsOklLokF8M3CAqxbu6Y1e0eYAXyDOtfACtj+fhu59GxVevZZ\nAgyASAIMgEgCDIBIAgyASAIMgEh5y+gL8wWAqvRsVXr2WUZgAEQSYABEEmAARBJgAEQSYABEEmAN\ncV+1qvRsVXr2WQIMgEgCDIBIAgyASO7E0Yy+951+gP2MwNowvQ7smjDADkZgDVgmVt937rFWiDvm\nVaVnn2UEBkAkAQZApLwpxNfXBpeD940fAVBMWID1ff8Kp+nrl/Gd1R+1axjml8GCKg/wkLAA25YU\nWjO5NQd4SKkAG82GX7M7lbUccmNVW64hn9GzVRXr2bj7OpYKsNWdqcy+BXCpjVP/NlVbhSiuAL5E\n2AhsGIblUsNxznB830JEgC8RFmDdWjKN70gsgK9SbQoRgC+RNwIrzCCyKj1blZ59lhEYAJEEGACR\nBBgAkQQYAJEEGACRBFhD+r6PuH0LR+nZqvTsswQYAJEEGACRBBgAkQQYAJHcSgqaMV0O4B5F8BMB\n1hD3VatqV8/OFrP1vQxrn8/ss0whQqusz4ZNAgyASAIMgEgCDFrl+gpsEmDQhllcSS/4iVWIDRlv\nqmZdUz17e1bXp/GZfVbZAHvdYfP14t1OtrwX58buuH/jQ8Xud1EF7in28QpcVOzNPXto48f766Ji\nH6/AIaG7d/vx3LdcuY/9P7TGf/59f9nY/VteV+zst/bsWO23a7blRcU+0l9/uun1z82N//mVI0ei\nHyuwp7YXtatqf33wsd1ziA/6cx362Lag1DWw/q+u64ZJH7xevz1M7Nhytdh31ThagT0S2zW8eX1W\nsfe3a3Y03Cj2FA22687+OrHY5W/9vgJ7ZO3ey6JaTK2JOgHW9/0wDIdOE5Y79+hdEe+2bO2BQI+3\na7sCFxW79MsK3Pzn2i+rXb/pr4uKvfljm9Kuiz62l6ozhTgG2J8Xqxu8+cX9G19U7H5Z7QqqwHXF\nXuHxdr092b+rAhcVe6hd+z3errOKbTAsyi7ieGfaB9vnFPu3/LjYQx6vbeMVuKjYjyuwX1a7Gi+2\n6sf28T9sm+pMIU6d3g/nrERqoITugtOoZa22/xMrl47PqMbj7bpIg+1a7a+jnXjRbrD9n9h2Sruu\nKGHmosPR0Y9tC2oG2E5/phwX77/b4d5t2VrXNtKu/RXY6aJ2vSv2ndP/XL8/zrbcrt/sBofaFfSx\nbbldp39sL1XnGlj30QzPZ8uXzy12v4sqcGexj1fgomJv69lDG+uvdxv/sgKHVP1ztaBUgE0tu2H/\nt3Cu2HJ74/2y2tVasTu/r3NusfsFtatqf93zsW28XRft3lcoG2BxpjtNmU55LQ3tJg0s0LrVdnXJ\nTVv2To3+2mhXl9y0qv111NetQmxZpT1vdhI3O+LntjRrguWQP9dg+n4YhjL91f3bruk76ar21yFf\nvYijNX3fJ65kXTUMQ8lP0Wq7CnRcyc7q3l/m0V81GIE1ZHmqSIQyHVegCatm7arRX+kZfAoB1oro\nz9I3q9Fx49GwRlumlu0q08ZXDD9dkSeZQmzCl++FuSp1XJkj+8y0XTX6q0YrTpE9iK6k5CKir1qF\nmNuu1WXT2tUyqxBHAgyASKYQAYgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgk\nwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTAAIgkwACIJMAAiCTA\nAIgkwACIJMAAiCTAAIgkwACIJMDgbn3fP10FqECAARDpv09XAL7RaxA2DMOzNYFcAgwe8Mqtvu9l\nGHzGFCIAkQQYAJFMIcIDXAOD3zP/DkAkU4gARBJgAERyDQyuNb3vhhl7OJFrYHCVd7eM8qGDU5hC\nBCCSAINLbNyx18184RQCDIBIxa+BOdUF+Ez76VB/FWL7ffBS9b6u39mu7ZOnlv8g39lfuS5qV8TZ\nvylEACIJMLjExklxyXEA3E+AwVVWg0p6wVlqTgq/VJ31Jst4OcGuSJCIg6cRGACRBBgAkfKW0W8/\nCdBcDcCXCAuw6bTscor29U7E7C0AvxEWYB+YfR1PsAGsivjy8lS1AFtOIUosHmcnJMJ0R40IszoB\ntj27CEAxViECEClsBDYMw3IV4jjeWv0RAFWFBVi3Fk6vd+QWwPcwhQhAJAEGl+v7PmJNF2QRYABE\nEmAARBJgAEQSYABEEmAARMr7HhjE8Q1FuIIRGACRBBgAkQQYAJEEGACRBBgAkfJWIW48M8XjVGjT\n8kHhwO+FBdjGY5enxwhPZAYoLyzAfrQ81Z3dBVywAayKe2ZCqQBbHZxJLIA9Nk7922QRBwCRSo3A\noDnjnPbkdWdKAE4SFmDDMCyXGo4Thqs/AqCqsADr1sLJ5S6AL+QaGFxm9TJ4wrVxiCDA4DKrUwLm\nCeAkAgyASHnXwCDJON6y/hAuYAQGl+u7zoUvOJ0AAyCSAAMgkgADIJIAAyCSAAMgkmX0cDk3OYMr\nGIEBEEmAARApL8D6vzY2uLM+ADwi7BrY+Oiv5evpBrdXCoAHhAXYtjHSZhk2+6fL6QCr4gYApQJs\nlcTiceNxwa5I46a7aESY5V0De2f8c0//H4DC6ozAtq+NAVBMWIBNL3FJLIBvFhZg3dqFhNk7wgzg\nG9S5BgbAV8kbgUEcswJwBSMwACIJMAAiCTAAIgkwACIJMAAiCTC43PYDgIDPCDAAIgkwACIJMAAi\nCTAAIgkwACLl3Qtx+TiVPT+CB9kh4QphATZ99NfqY8DGdzwhDKC8sADbthpas+/fCDaAVXHfViwV\nYKPZ8EtiAewxPVpGhFmpABv/4hIL4BtUW4UovQC+RNgIbBiG5VLDcc5wfN9CRBpkbgCuEBZg3dpR\nYHzH0QHgq1SbQgTgSwgwACIJMAAiCTAAIuUt4oA4VhjBFYzAAIgkwACIJMAAiCTAAIgkwACIJMDg\ncn3fRzycArKUXUb/Ol78eHvf5ZFlY9Hz/o0PFbvfRRW4p9jHK3B/sftltSuoAtcVu19Wu2YHz5a/\nBNK3XLmP/fm7v/7598Wysf+Pt3833ni48+nF7hfUrndbbmxcr13zmhw5Em2U+Xi7Zht/UIGNjSPa\ndUj53ftBpaYQ+7+6SQeMr4e/G8y2n/70tfFyy26xu8xef1zsfvvbtVrDQ+36/Z/rXQWWqrbrFA22\na1nsB/21Uex0y8bbtUfi4ei23fv36gTY+FSwU04T3hWxfP/3+/cVljv3aH+7Lq3Ax3+um9t1225w\n9Bez2lWvv45KaddFH9tL1ZlCHAPsz4vVDd784v6NLyp2v6x2BVXgumKv8Hi73p7s31WBi4o91K79\nHm/XWcU2GBZlF3G8M+2D7XOKWW/t3/hQscv/xP7LsKdU4IpifzxZu6LY1tq1X1a7xjn6ZTuXx7jH\nd4OL2vXPTz+9unm0ti3sBg2qM4U49ct+uKgbf19smxVb/vr2cXzl0vHvKvDOze26TY12nb4bNNKu\npZTD0dGPbQtqBtjS6mnUnynHg4WcXuxvHKrA0ZZe9Oc6vQKNtOvjT/gJS/DXymmkXfX66/daPhz9\npr/uV+caWPfRDM8Hy5dPL3b2W6cstm6hXasbP16Bi4o9d6V1t3mQfbxdjxf7eAVmv7Lnc72n2Kq7\n93VKBdjUrBsO7WH7I+SsYqfbn1jm4+0KqsB1xe61PHC03a7pxo9X4KJiT/nM/qYC2xu3UOyzygbY\n6LU0McLODwM1rZ4m2xnaVvgzG3Hw/LpViG2anvJE7DdAVzS6gnzLIo6WvfuaPQAbBNjDpjd6ef2v\nk2F05g/hBwLsecujlOPWN5rFlfSCn7gG9qTtYZaLYV9Hd8MRRmAARBJgAB/q+97l6gcJsCdt32vH\n/CHABgHWhP7NawDesYjjYX8GYf9OQxh7AfxIgDVBYgEcZQoRgEh5I7DXZNvGUwMMaIAbONQ8KyzA\npt/tXf2e7+uSkh0LoLawANu25/lsgg1gVdx32koF2Gg2/JJYAHtMj5YRYdZ0gF3xRGMAamg6wD6I\nIukF8CWaDrClYRiWSw3HOcPxfQsRgduY9XlWWIB1a/vK+I59COCr+CIzAJEEGACRBBgAkQQYAJHy\nFnEANMLasWcZgQEQSYABEEmAARBJgAEQSYABEEmAAXyo7/uIx45UJcAAiCTAAIiUF2D9Xxsb3Fkf\nAB4RdieO8dFfy9fTDW6vFPC93IrjQWEBtm36ZMvpm9N/uvULcI6+f51Nd13X5R9b4gYApQJslcQC\nzrc81vd9eoZNj5YRYdZ0gB0aPI0bv/5fbgHU1nSAHQqh7WtjABTTdIAtTS9xSSyAbxYWYN3asGz2\njjADLjcM88tgjjy3y/seGEATXok1DNLrEQKsIRHLfj6gXVm0a2dp4/0U+qfviFi1v/YQYADHrGaG\nG/veT4ABEEmAARxgmNWO4gvQ7WoAn2k/HYoHGMC5fjwtdlC9jSlEgNNIrzsJMIADRFQ7BBjAMe8y\nTLbdTIABHPbjPe24gUUcrZheGS7TKbMnaI8vCrRutV1dctNW75E9eyfRRru65KZV7a+j8m7mW1il\nPW/5XOwaTw9YXYGW25ypsRWvx5rX6K/u33ZN30lXtb8OMYXYkEq3ohmGoeSnaLVdBTquZGd1b9ql\nv8owAmvI8lSRCGU6rkATVs3aVaO/0jP4FAKsFdGfpW9Wo+PGo2GNtkwt21Wmja8YfroiTzKF2IQv\n3wtzVeq4Mkf2mWm7avRXjVacInsQXUnJRURftQoxt12zA+Ls1F67GmQV4kiAARDJFCIAkQQYAJEE\nGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgAkQQY\nAJEEGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgAkQQYAJEEGACRBBgA\nkQQYAJEEGACRBBgAkf4HaNN6Ul9C7KIAAAAASUVORK5CYII=\n", "output_type": "display_data"}], "prompt_number": 28, "cell_type": "code", "language": "python", "metadata": {}, "input": ["subplot(2,1,1);\n", "hold on;\n", "mystem(1:N, a0, 'k');\n", "plot(a, 'r.', 'MarkerSize', 20); \n", "axis tight; box on; title('a');\n", "subplot(2,1,2);\n", "hold on;\n", "mystem(1:N, b0, 'k');\n", "plot(b, 'r.', 'MarkerSize', 20); \n", "axis tight; box on; title('b');"]}, {"source": ["__Exercise 3__\n", "\n", "Compute the full homotopy path $\\la \\mapsto (a_\\la,b_\\la) $\n", "where $ (a_\\la,b_\\la) $ is the solution of C-BP with regularization parameter\n", "$\\la$. The display bellow shows in color the evolution of the correctly estimated spikes\n", "and in black wrong perturbating spikes.\n", "Test with different values of $\\kappa$. What can you conclude about the performances of C-BP ?\n", "ist of regul param\n", "armup stage\n", "uter loop\n", "isplay"], "metadata": {}, "cell_type": "markdown"}, {"collapsed": false, "outputs": [{"metadata": {}, "png": "iVBORw0KGgoAAAANSUhEUgAAAkAAAAGwCAIAAADOgk3lAAAACXBIWXMAAAsSAAALEgHS3X78AAAA\nIXRFWHRTb2Z0d2FyZQBBcnRpZmV4IEdob3N0c2NyaXB0IDguNTRTRzzSAAAYtElEQVR4nO3dbXfa\nuBqGUZjV/z1zfjnnAy11eDGWZVu65b0/zOokaSKs9LliQ+B6u90uAJDmn9YLAIA1BAyASAIGQCQB\nAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBI\nAgZApF+tFwBndL1eH3++3W4NVwK5BAwaeERrWjKgiIBBA7oF9QQMjna9Xp2BQT0BgwZ0C+pd3YEM\nQCIPowcgkoABEEnAAIgkYABEEjAAInkYfYUVj4T2mE+AjQwesP1+22ZliK7Xy+XiN4CA/vX/S1aD\nB+yyxx5Mo7jgk/9+2oXpk7eW/PV60+d9aKL5AnpYgwX0sAYLWL6GiN+1dx9YocJ6/fjg14+/Xtdc\nhwTgDGdgW1pdr9e/Nf1U9z+3/rkMIIszsMXq6zX1ekL282zs+t/1+p+TM4CP2l+Q3dVmV5y3rdfM\nJ7+/4d+PH3v7b+T9AjrRw911XwUsscY2e7Brvd5+lcfbPpfsQdKAzQlYexvswSH1elwtvP3v3btv\nt6LLiZIGVBKw9qr24OmUaJ8D9Vqm3/l5PSH7uYDlSdMzoJSAtbd+D/Y/8fqYrk/LmF2MngEbErD2\nVu7BzvValK5P6/n9F76samHPxAx4S8DaW7MHe9arOF0//vK7Ji1boZ4BRQSsvev1Onfz3j41xsx7\na1ZSk64fn2h9xmYW85aewWkJWHtfAjZjz3ptEIYtMvb3k+kZ8JOAtbcyYNsdk+3T9fwFiu8e+/4p\n9QxOT8B28XiO5NeVv76r7R7sXq+/X2n7jP393HoG5yNg25se06fje6/X/S2Pd7Xag+PS9eOr7pix\nv19Ez+AEBGx7XwN2N/2Yp89wwO1tU6+/X37Lu8e+fzU9g1E8Dcz+6zBUwF7fdfAPEY3TNXVsxv5+\nWT2DITgD217PAeuoXlOHXFf8+MX1DDIJ2PZWBOxyuV12Htqdpmuqacb+rkLPIISA7WLmoYZv33UP\n2NS2tzigXg99ZOxBz6BbAtbe23uCpmpufVK6pjrL2J0nI4auCFh7P68rfvngoiORWq+HRo/yWMjJ\nGbQlYO192oP5mM0fkvh0TfWdsQc9g4MJWHtL9uBTzN480+9I6XrS5XXFT/QM9iZg7ZXuwacxvtlz\nyXcuKmMPegabE7D21u3BjzH+7znSNZWZsQc9g3oC1l7NHryZg/8b+2j9FJ6xBz2DFQSsvZVnYO/S\nNf2/oY/ZTyGP8lhOz2AJAWvver1e/q36DPdBNtwYLzTu7dczeEvA2qsJ2NuZNcp1tbVGv/1LeiZm\nnIGAtbcuYF8n1Ohj/Jtz3H4nZ5yZgLW36x6cY4x/drLbr2ecioC1d8AejHv30DJnvf0uNjI2AWvv\nsD046xj/4+y3X88YjYDt4vU1U17f9Xjv8XtwsutqL85++y8XMWMIAra9mRe0XPiWY5x9jJ/99v/w\ntWdiRocEbHtfA3b/w8zHHOnsY/zst/8NMSOFgG1v4RnY62s0Pxx/e88+xk9/99gnYkZvngZm/3UY\nKmCvH9bPDxFnH+Nnv/1fzMdMyTheP8NzRsASp2YC9vZdve2BMX76E9LvZmKmZBymt+H5VsASn8zc\n0dXbfWAzzj7GlXwZMaOVbofnVMASa3S+BzL25o3nOgRLKRkH63x43gUssUbGHsjYW+c6CkspGcfI\nGJ79L7FGxB7cORtxCEq9jZmMsYmI4RmwxBoRezBlhjsEKygZm4sYngFLrBGxB2+d/briRcmKyRgb\nihieAUusEbEHM2RMxlZQMupFDM+AJdaI2IOvZOxyUbJiMkaNiOEZsMQaEXuwkAF+uTgKaygZK0QM\nz4Al1ojYgyIG+G8ORCEZo0jE8AxYYo2IPVjB9P7NgSinZCwRMTwDllgjYg9quHvsNyUrJGPMixie\nAUusEbEH9WTsNxkrp2S8FTE8A5ZYI2IPtiJjfylZIRnjScTwDFhijYg92JbR/ZdjUU7JuIsYngFL\nrBGxB3swun9wOArJGBHDM2CJT15f9Ov1A/p/PbDDuK74l4wVkrEzixieAUuc+hqne94E7ImM/aBk\nhZTshCKGZ8ASp+YDdn+LM7BPZOwHGSskY6cSMTwDlji1ImBPnyHr9u5Bxp4pWSElG9XTwOx/Wo4T\nsLeHPuKHiCYM7WeOSCEZG1vE8AxY4tSSy4MuIS5naD97e0QuJz8oXyjZkCKGZ8ASn7w+CvH1VEzA\nSrmu+EzbC8nYYCKGZ8ASa0TsQT9k7JmMlVOyMUQMz4Al1ojYg97I2BtKVkjG0kUMz4Al1ojYgz7J\n2BsyVk7JQkUMz4Al1ojYg56Z2O85LoVkLE7E8AxYYo2IPeifcf2e41JOyVJEDM+AJdaI2IMgriu+\np2SFZKx/EcMzYIk1IvYgjoy9J2OF3mbsomR9iBieAUusEbEHoWTsIyUr5ISsQxHDM2CJNSL2IJpZ\n/ZFDU8gJWVcihmfAEmtE7MEAzOo5jk4hJ2Q9iBieAUusEbEHI3Fd8SMZK+SErK2I4RmwxBoRezAe\nGZujZIWckDURMTwDllgjYg9GJWNzZKyQjB0sYngGLLFGxB6MTca+ULJCSnaMiOEZsMQaEXtwBqb0\nFw5QIRnbW8TwDFjik9fXA5t5V8QenIcp/Z1jVEjJdhIxPAOWODX/YpWPt7z+ga64rviFjBWSsc1F\nDM+AJU4tfLXlacCe3pV1e8cmY1+8zdjFYZqjZDWeBmb/03K0gD1dRXwN2JOsmz8kGfvOCVkhGavn\nDGx7K87Abrfb14xNZR2QYcjYdzJWTslWE7DtzQTs7bs+7cGSpGUdmTEY0Ys4TIVkbAUB28XMQw1r\nHoU4n7S4oxTNfF7EYSqnZMsJWHvr9mAmZmMfrt64rriIkhWSsSUErL3KPVCyHsjYIjJWTslmCFh7\nW+2BkjUnY0spWSEZe0vA2ttjD/xuWUOG81KOVDklmxKw9nbdg7hf+huG4VzAwSokY3cC1t7ee6Bh\nbbmuuJSMlTt5yQSsvWP2YJqxsY9nn2SsgJIVOm3GBKy9w/bAqVhzMlZAxsqdrWQC1t7Be+BUrDmT\nuYzjVeg8GROw9o7fA6diPTCWyzhehd5m7DJWyQSsvSZ7oGH9cF2xjJIVGviETMDaa7gHLif2Q8bK\nyFihITMmYO213QOnYl2RsWJKVmikkglYez3sgVOxrshYMRkrNEbGehieXwUs8cnra6bMvKuTPXAq\n1hszeQ1HrVB0yToZnvMCljg1/4rMr69j2c8eaFiHDOQ1HLVCoRnrZ3jOCFji1HzAXj+stz1wObFP\nrisW+/T6DA7cZ1kl6214vhWwxKklAXv6mKf3Nr+9TsW6JWNrOCEr1HPG4qbTUAG7H/3pG7v9IcKp\nWLdkbA0ZK9dzyS4dD8+pgCVOLbwPbOYt/Yj7YedUZGwlJSvUbcZ6Hp4PAUt88umhhm970PkeaFjn\nTOOVHLhyvZWs8+F5F7DEGhl74HJi30zj9Ry7Qv1kLGN49r/EGhF7cHEqFsJ1xZVkrFAPTxYcMTwD\nllgjYg/uNCyFjK2nZIUanpBFDM+AJdaI2IMpGUshY+vJWKEmGYsYngFLrBGxB09ef3ftLu6GnIFR\nXMXhK3RkySKGZ8ASa0TswVufMvYQeruGZA5XcfgKHZOxiOEZsMQaEXsw72vJLmLWDdcVqyhZoV1L\nFjE8A5ZYI2IPiuhZ/2SsiowV2iljEcMzYIk1IvZgtSUxu+hZIzJWxZMFl9u2ZBHDM2CJNSL2YCt6\n1iEZq+WErNBWGYsYngFLrBGxBzvRs34YwrWckJWrLFnE8AxYYo2IPTjAwphd9GxPMrYBB7HQ6oxF\nDM+AJdaI2IPj6VlbrivWckJWrrRkEcMzYIk1IvagORcbm5CxDTghK7Q8YxHDM2CJNSL2oDd6diQZ\n24ATsnJfSxYxPAOW+OT19cBeP2DmRS8ppWcHcCKxDcex0EzGIoZnwBKnvr4i8/0PArYfPduP8bsN\nx7FcP69DVuRX6wVs6fEqzNM33v/XPN3Kpx8annha/RXuB+npiN7/1/Er8DhY00P5+LND+c7vs64P\nr0PWraECNsOrHu9EzzZn/G7GTwSl/vfnD/+2XMVyYVfYlty/9fQx858w6+bHcb2xnuthm3EoS0Tc\n/3KWM7BPvj4khBrOz+o5IduME7LhBDT2yddHasycpTkh6IrtWMFZxGYcylkRZ2ABS6wxsweuLvZG\nz4oYv9twHD8QsPYW7sHM6Bz7+PRMz5YwfjfjUP4kYO2VPvjw09Ac+yhF0LN5xu82HMc/BKy9+an3\n6UGMRR9PE3r2lvG7mdMfSgFrb/nTrl9+PipkyYfRDz17cvrxu5ETH0cBa+/18dlLJt3bZ/T49GH0\nxovFPJx4/G7tfIdSwNpbci41c7/X8trRLT27nHH87uNMx1HA2rvvwfLrS1/P2GaaN/aRHMaZe3am\n8bunc7x6i4C193YPiu4Yu7xkbP6usrGP53hO2LNzjN9DDP0TgYC1t2QPFl4n/PSIfI+8H8mpejb0\n+D3QoD8RCFh7pXuw7lGL83937CM8vOEf3zjo+G1hrJ8IBKy91XtQM7ZcWhzYwKdoY43fdkb5iUDA\n2qvcg8pppWTDG69no4zfDoT/RCBg7W21BzVzqvlTezT/Rmy+gMPWUHMJem/ll9Pfv71m1c2/Exos\noL+fCJY/Q2zzf7Zfjf96YPeZUrkTC3+1efox06/46a97NbLBFP348vqurr4NHmt5+/pZl6RziaYc\nxz0NFbCZHmwyLJa/0PP8V5wp2bqF0a2nl6ab/87p83FAJvA2/vwU+zqbHu+l1DgBm3kdy08f//SW\n5WNi+QnZ9Gs9Pv/8o/D7nGJsYt1vJTa/Cv3na/3+g5JVeXscHcRVAq5yLvQ2YKW/swzAXf91+Kf1\nAvbV/wYAsM44lxA/0TCAIQ1+BgbAqMa5D+ziUekAZzJUwAA4j/HvA2vo9cGQh/248PoVD15D8wVM\nV9J8Fy4vv3pxkl2YeYG9I78Nnn6P5YS78OnbL/2qlYDtYvrvtvQX1Lby+Ja9T43j19B8AZfJRrTd\nhSZrmE7tVrvwdPOb7MJ0fB9/EF6/3JELmB9Erf5RbMiDOHZxu91aP3tC4+/F5gu49PFv8j61my+g\nh+PQfA0n1HwQ7c0Z2MjaTg2/RX75eQLU6qs3XEAnNnlC1E3WwLYEbEw9/Itt+2Qo9687/e/xztyM\nqYb5bH6V7HG339dnwmQFlxCHdfJzr9sfl0aHooeDQHP3avpRZienvrCwt1aPf3t96FfDNTR/vFMP\nj0I87ePfLi/nPXahyQJGfRSigAEQySVEACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAi\nCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAA\nRBIwACIJGACRBAyASAIGQCQBAyCSgMHRrtdr6yXACAQMgEi/Wi8AzuhxEna73dquBHIJGDTw6Nb1\netUwWMclRAAiCRgAkVxChAbcBwb1XH8HIJJLiABEEjAAIgkYAJEEDIBIAwbME80BnMFQD6OXLoDz\nGOoM7Ha7+a0AgJMY6gzslXMygHX6Px8YPGCXDvag+bO1WkAPa7CAHtZgAcvXEPHT/1CXEAE4DwED\nINKAAWt+hg7AAdpfkN1VD1ecAeJEDM8Bz8AAOAMBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgC\nBkAkAQMgkoABEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDINKv1gvgjet/1+n/3v7r/YW9AY53\nvd1GHo7Xa8ANfMrVJjQPqJExPPtfYo0+92CPYq0mdcCrPofnk4Al1uhkD5YXqzQnR7ZQ6uA8Ohme\n8wKWWKPVHizsypFJOPi0T+0gmoC1d+QefC1E0ExXOzg5AWtv7z0YKVrrSB0MScDa23YPOrwwmMU9\ndpBCwNqr2YOiaWtcbsuJHbQlYO193YPVg9LI64cTO9icgO3iev09rV5X/njX473X6/Xy7wZf1Nga\nw2Gp8w1DOgHb3vSYvh7ft29ZETDTh2NS5zuNbgnY9r4G7P6HmY+Bbe2aOoWjlYjhGbDEqYVnYNM/\nPH2GrNvLMDzjJf17Gpj9T8uhAvb6YRE/RMDdVpETNupFDM+AJU7NBOztuyL2ABaqLJywsVzE8AxY\n4pOZO7rcB8Zp+YUQthUxPAOWWCNiD2AnqsZqEcMzYIk1IvYADrYubKp2KhHDM2CJNSL2AHqgakxF\nDM+AJdaI2APolouQpxUxPAOWWCNiDyCLqp1BxPAMWGKNiD2AMbgIOZKI4RmwxBoRewADU7VQEcMz\nYIk1IvYAzsZFyP5FDM+AJdaI2APgomqdiRieAUusEbEHwCc1z54lbDUihmfAEmtE7AFQyuna3iKG\nZ8ASa0TsAbAJVdtQxPAMWGKNiD0AduWRkCtEDM+AJdaI2APgeE7X5kUMz4Al1ojYA6ATqvYQMTwD\nllgjYg+Anp2zahHDM2CJNSL2AEg09l1rEcMzYIk1IvYAGMYwVYsYngFLrBGxB8DYEi9CRgzPgCXW\niNgD4IQ6r1rE8AxYYo2IPQB46OQiZMTwDFhijYg9AJh3fNUihmfAEmtE7AHACrtehIwYngFLrBGx\nBwBb2apqEcMzYIk1IvYAYG8rwtbhg/uf/Gq9AAB297ZGNS+31oPxA3a9VrwgnrM3YFzpVRv8CltN\nvd4a+3AB3EXc/zL+Gdi23hax/20GGI+AbeDTeZ6wAexnqIA9QrK8HK8fueFVx5lPpW0AlcYJ2PSK\n7cKrt28/5umNm9+L9vXTahvAEuMEbIWd4lSpz1UB9Oaf1gsAgDVOfQZ2+XC9bo9zIBcG0/VzZrzH\n91L9jfMN/lvNoezpIPbzDT9j/IDdbreZndh2k1QqzvH/Spt/k4wyYDuw+lA6jhsZJ2DTUE1nxE4T\nqvkY4q0TBmmGAbsNx7FX4wTsUjFKpg9f3OPzU+ngLMVttAG7GYcyylABW82D2ps4MkvD7KMBuxmH\nMp+AvTfMvDuYU6WtmK5bWnc0HcrunSFga74LEx6Ac5jjj8XSLTv5Nhmwz2T/ZM4QMN46ePYbEFUM\n2GdaxRkCdqpvV/cqMSCt4oPxAzaGw8okSzQjVBQSsL74rTVOQavYgoC1tzpaskQArWI3AtbG12iJ\nE3m0imMJ2HFEi0EIFX0QsH2JFtm0io4J2PZEi0ieroI0ArYN0SKJVjEEAVtpyUMHRYvGXABkaOMH\n7FNpltSl9AHuikUzWsX5jB+wT7b6lWHR4mhaBZfL5cwBW0GrOJo7q+Cz8QOmOgRwUgXlxg8Y9EWr\nYCMCBrvRKtiTgMEW3FkFhxMwKOGkCrohYPCBVkHfBAy0CiIJGGciVDAQAWNQWgWjEzDyeQQgnJKA\nEUWrgD8EjC65AAh8kxewx7PIvz7J4fQJ5j0FYgytAlYJC9j1en2Uafrnhy66tdELtfBbD3sK9Ccs\nYF/dT8J+ZGySk+u6wgnSYbQK2tnqVRIPM1rA7n16e3J2uVzu7zt6TX+/vOkM9Gs6NiNi1nXAno7g\n15OnBtcPNQmgka4DVhSkT2ddGgMwpK4D9up2u70+CvGerrfvAmBUYQG7vIvT4y26BXAe/7ReAACs\nIWAARBIwACIJGACRBAyASAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoAB\nEEnAAIgkYABEEjAAIgkYAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJ2O6u16sFtF1AD2uw\ngB7WYAGdrGErAgZApAEDNtLPFwB88qv1ArYkXQDncb3dbq3XsLHr9e+NkjSAdfqvw1BnYK/63wAA\n1kkN2NOplVABnE1qwBQL4OQGfBQiAGcw4IM4ADiD1EuISzzuJ2sV6cfjIY9fyetXPHgNzRcwXUnz\nXXh80bPtwut91U2+De5ftPlBaLiAT99+zYdkpWED9vRg+uOHZvOVPL5l71Pj+DU0X8BlshFtd6HJ\nGqZTu9UuPN38JrswHd/HH4TXL3fkAuYHUdshuQn3ge3idru1/W5o/r3YfAGXPv5N3qd28wX0cBya\nr+GEmg+ivQ17Bsal9dTwW+SXnydArb56wwV04ukSYsM1sC0BG1MP/2KfLrsf7P51p/893pmbMdUw\nn82vkj3u9pve/8dWXEIc1snPvW5/XBodih4OAs3dq+lHmZ2MfGGh+QNsWj3+7e3TlJztkVfTlTR/\nFOJpH/92eTnvsQtNFjDqoxBHDhgAA3MJEYBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyA\nSAIGQCQBAyCSgAEQScAAiCRgAEQSMAAiCRgAkQQMgEgCBkAkAQMgkoABEEnAAIgkYABEEjAAIgkY\nAJEEDIBIAgZAJAEDIJKAARBJwACIJGAARBIwACIJGACRBAyASAIGQKT/Ax+ZZfKOPcMIAAAAAElF\nTkSuQmCC\n", "output_type": "display_data"}], "prompt_number": 29, "cell_type": "code", "language": "python", "metadata": {}, "input": ["exo3()"]}, {"collapsed": false, "outputs": [], "prompt_number": 30, "cell_type": "code", "language": "python", "metadata": {}, "input": ["%% Insert your code here."]}, {"source": ["Bibliography\n", "------------\n", "<html><a name=\"biblio\"></a></html>\n", "\n", "\n", "* [Ekanadham11] C. Ekanadham, D. Tranchina, and E.P. Simoncelli. _Recovery of sparse translation-invariant signals with continuous basis pursuit_. IEEE Transactions on Signal Processing, 59(10):4735-4744, 2011.\n", "* [CombettesPesquet10] P.L. Combettes and J-C. Pesquet, _Proximal Splitting Methods in Signal Processing_, in: Fixed-Point Algorithms for Inverse Problems in Science and Engineering, New York: Springer-Verlag, 2010.\n", "* [DuvalPeyre13] V. Duval and G. Peyre, _Exact Support Recovery for Sparse Spikes Deconvolution_, preprint hal-00839635, 2013\n", "* [deCastroGamboa12] Y. de Castro and F. Gamboa. _Exact reconstruction using beurling minimal extrapolation_. Journal of Mathematical Analysis and Applications, 395(1):336-354, 2012.\n", "* [CandesFernandezGranda13]  E. J. Candes and C. Fernandez-Granda. _Towards a mathematical theory of super-resolution_. Communications on Pure and Applied Mathematics. To appear., 2013.\n", "* [Chen99] S.S. Chen, D.L. Donoho, and M.A. Saunders. _Atomic decomposition by basis pursuit_. SIAM journal on scientific computing, 20(1):33-61, 1999.\n", "* [Yu13] Y. Yu, _On the decomposition of the proximal map_, Proc. NIPS, 2013."], "metadata": {}, "cell_type": "markdown"}], "metadata": {}}], "nbformat": 3, "metadata": {"kernelspec": {"name": "matlab_kernel", "language": "matlab", "display_name": "Matlab"}, "language_info": {"mimetype": "text/x-matlab", "name": "matlab", "file_extension": ".m", "help_links": [{"url": "https://github.com/calysto/metakernel/blob/master/metakernel/magics/README.md", "text": "MetaKernel Magics"}]}}}