{"id":10035,"date":"2018-09-07T09:00:05","date_gmt":"2018-09-07T13:00:05","guid":{"rendered":"https:\/\/blogs.mathworks.com\/pick\/?p=10035"},"modified":"2018-09-07T13:35:18","modified_gmt":"2018-09-07T17:35:18","slug":"path-simplification-and-binary-image-reconstruction-made-easy","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/pick\/2018\/09\/07\/path-simplification-and-binary-image-reconstruction-made-easy\/","title":{"rendered":"Path Simplification (and Binary Image Reconstruction!), Made Easy"},"content":{"rendered":"<div class=\"content\"><!--introduction--><!--\/introduction--><\/p>\n<h3>Contents<\/h3>\n<div>\n<ul>\n<li><a href=\"#a88e24bd-6931-4760-b2c0-1a8bda5472e5\">Sample Path<\/a><\/li>\n<li><a href=\"#f96d4210-5955-43cc-a179-714615fb820c\">Enter the Douglas-Peucker algorithm:<\/a><\/li>\n<li><a href=\"#93121bd8-d8d0-481a-a9dd-861edc2bd4be\">Mask reconstruction<\/a><\/li>\n<\/ul>\n<\/div>\n<p>I&#8217;ve been thinking a lot about how to generate and efficiently store paths&#8211;and then reconstruct binary images from them. Path simplification is important in many fields, and after some Googling, I found my way to the <a href=\"https:\/\/en.wikipedia.org\/wiki\/Ramer%E2%80%93Douglas%E2%80%93Peucker_algorithm\">&#8220;Douglas-Peucker&#8221; algorithm<\/a>. (It goes by other names as well, but that seems to be the most common.)<\/p>\n<p>It turns out that Douglas-Peucker is used extensively in cartography to simplify the representation of geographical borders. Armed with that information, I came back to MATLAB&#8217;s doc-search engine and, by typing &#8220;Douglas Peucker,&#8221; I found my way to the <a href=\"https:\/\/www.mathworks.com\/products\/mapping.html\">Mapping Toolbox<\/a> implementation of <tt>reducem<\/tt>. From the doc:<\/p>\n<p><i>The <a href=\"https:\/\/www.mathworks.com\/help\/map\/ref\/reducem.html\"><tt>reducem<\/tt><\/a> function implements a powerful line simplification algorithm (known as Douglas-Peucker) that intelligently selects and deletes visually redundant points.<\/i><\/p>\n<p>For those of you (us) who have the Mapping Toolbox, you&#8217;re good to go. Despite the fact that it was implemented for a non-Cartesian geometry, it works nicely, as I will demonstrate. For others, a quick check of the File Exchange reveals <i>four<\/i> implementations of this particular algorithm:<\/p>\n<p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/pick\/files\/FexSearch.png\" alt=\"\"> <\/p>\n<p>Let&#8217;s create a sample path, and consider each of them (along with <tt>reducem<\/tt>) in turn:<\/p>\n<h4>Sample Path<a name=\"a88e24bd-6931-4760-b2c0-1a8bda5472e5\"><\/a><\/h4>\n<pre>pgon = nsidedpoly(5,'Center',[150 150],'SideLength',100);\r\nverts = pgon.Vertices;\r\nverts = [verts;verts(1,:)];\r\nverts = [verts(1:2:end,:);verts(2:2:end,:)];\r\nmask = poly2mask(verts(:,2),verts(:,1),300,300);\r\nmask = imclose(mask,strel('disk',4));\r\nmask = imfill(mask,'holes');\r\nfigure\r\nsubplot(1,2,1)\r\nimshow(mask)\r\ntitle('BW Star')\r\nperim = bwboundaries(mask);\r\nperim = perim{1};\r\nsubplot(1,2,2)\r\nplot(perim(:,2),perim(:,1),'b-')\r\ntitle('Extracted Boundary')\r\naxis equal\r\nset(gca,'xlim',[0,300]+0.5, 'ylim',[0,300]+0.5)<\/pre>\n<p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/pick\/files\/star.png\" alt=\"\"> <\/p>\n<p>If you run that code, you&#8217;ll see that variable <i>perim<\/i> has 547 [row, col] coordinate pairs. <b>Most of them are redundant for our purposes<\/b>, and we would like to be able to discard them.<\/p>\n<h4>Enter the Douglas-Peucker algorithm:<a name=\"f96d4210-5955-43cc-a179-714615fb820c\"><\/a><\/h4>\n<div>\n<ul>\n<li><a href=\"https:\/\/www.mathworks.com\/matlabcentral\/fileexchange\/61046-douglas-peucker-algorithm\">Douglas-Peucker Algorithm<\/a>, by <a href=\"https:\/\/www.mathworks.com\/matlabcentral\/profile\/authors\/2664900-reza-ahmadzadeh\">Reza Ahmadzadeh<\/a>:<\/li>\n<\/ul>\n<\/div>\n<pre>x = perim(:,2);\r\ny = perim(:,1);\r\nsimplified = DouglasPeucker([x';y'],0)<\/pre>\n<p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/pick\/files\/DPReza.png\" alt=\"\"> <\/p>\n<p>Wait&#8230;what&#8217;s going on there? Turns out that Reza&#8217;s implementation considers the distance between the first point and the last. That isn&#8217;t very useful for a closed curve! I&#8217;m going to move on to&#8230;<\/p>\n<div>\n<ul>\n<li><a href=\"https:\/\/www.mathworks.com\/matlabcentral\/fileexchange\/36229-douglas-peucker-function\">Douglas-Peucker Function<\/a>, by <a href=\"https:\/\/www.mathworks.com\/matlabcentral\/profile\/authors\/3289512-ranadeer\">Ranadeer<\/a>:<\/li>\n<\/ul>\n<\/div>\n<pre>[list_pts,pnts]= douglas_peucker([x,y], 0);\r\nsize(pnts)\r\nans =\r\n     2     2<\/pre>\n<p>Here again, the results match Reza&#8217;s (the two returned points are the point on the &#8220;upper left&#8221; prong, and are coincident). Ranadeer&#8217;s implementation also relies on distances from the endpoints of the curve. (Note that changing the tolerance in either Reza&#8217;s or Ranadeer&#8217;s implementation has no effect.) Let&#8217;s try&#8230;<\/p>\n<div>\n<ul>\n<li><a href=\"https:\/\/www.mathworks.com\/matlabcentral\/fileexchange\/41986-ramer-douglas-peucker-algorithm-demo\">Ramer-Douglas-Peucker<\/a>, by <a href=\"https:\/\/www.mathworks.com\/matlabcentral\/profile\/authors\/4185006-ligong-han\">Ligong Han<\/a>:<\/li>\n<\/ul>\n<\/div>\n<pre>simplified = DouglasPeucker([x,y],10);<\/pre>\n<p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/pick\/files\/DPHan.png\" alt=\"\"> <\/p>\n<p>Beautiful! That&#8217;s what I&#8217;m talking about! Conveniently, Ligong coded to allow an additional input flag to show the data:<\/p>\n<p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/pick\/files\/DPHan2.png\" alt=\"\"> <\/p>\n<p>Nice!<\/p>\n<p>Next up, <a href=\"https:\/\/www.mathworks.com\/matlabcentral\/fileexchange\/21132-line-simplification\">Line Simplification<\/a>, by <a href=\"https:\/\/www.mathworks.com\/matlabcentral\/profile\/authors\/870595-wolfgang-schwanghart\">Wolfgang Schwanghart<\/a>:<\/p>\n<pre>[ps,ix] = dpsimplify([x,y],10);<\/pre>\n<p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/pick\/files\/DPSchwanghart.png\" alt=\"\"> <\/p>\n<p>Again, that&#8217;s lovely! (Note that Wolfgang also outputs idx&#8211;an index of the points that are retained. That could be useful!)<\/p>\n<p>Finally, consider the Mapping Toolbox implementation. Even though it is implemented for geospatial (latitude\/longitude) data, it works beautifully!<\/p>\n<div>\n<ul>\n<li><tt>reducem<\/tt>:<\/li>\n<\/ul>\n<\/div>\n<pre>[xReduced,yReduced] = reducem(perim(:,2),perim(:,1),10);\r\nplot(xReduced,yReduced,'r.','Markersize',14)<\/pre>\n<p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/pick\/files\/reducem.png\" alt=\"\"> <\/p>\n<h4>Mask reconstruction<a name=\"93121bd8-d8d0-481a-a9dd-861edc2bd4be\"><\/a><\/h4>\n<p>So clearly we have some great options, even for simplifying closed curves. Once you&#8217;ve got that reduced set of vertices, note that you can readily reconstruct the polygon using poly2mask():<\/p>\n<pre><p>mask2 = poly2mask(ps(:,1),ps(:,2), 300, 300);<\/p><\/pre>\n<p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/pick\/files\/Reconstructed.png\" alt=\"\"> <\/p>\n<p>It is important to note that, while those binary images are nearly identical, they are not <i>exactly<\/i> identical. Why? Steve (in one of his stellar <a href=\"https:\/\/blogs.mathworks.com\/steve\/2014\/03\/27\/comparing-the-geometries-of-bwboundaries-and-poly2mask\/\">blog posts<\/a> on image processing) discusses that issue. It turns out that we <i>can<\/i> get identical results, if we upsample the path, with nearest-neighbor interpolation, by a factor of 3!<\/p>\n<p>As always, I welcome your <a href=\"http:\/\/blogs.mathworks.com\/pick\/?p=10035#respond\">thoughts and comments<\/a>.<\/p>\n<p><script language=\"JavaScript\"> <!-- \n    function grabCode_c0f6d247caae4ba8a0869999e3f60795() {\n        \/\/ Remember the title so we can use it in the new page\n        title = document.title;\n\n        \/\/ Break up these strings so that their presence\n        \/\/ in the Javascript doesn't mess up the search for\n        \/\/ the MATLAB code.\n        t1='c0f6d247caae4ba8a0869999e3f60795 ' + '##### ' + 'SOURCE BEGIN' + ' #####';\n        t2='##### ' + 'SOURCE END' + ' #####' + ' c0f6d247caae4ba8a0869999e3f60795';\n    \n        b=document.getElementsByTagName('body')[0];\n        i1=b.innerHTML.indexOf(t1)+t1.length;\n        i2=b.innerHTML.indexOf(t2);\n \n        code_string = b.innerHTML.substring(i1, i2);\n        code_string = code_string.replace(\/REPLACE_WITH_DASH_DASH\/g,'--');\n\n        \/\/ Use \/x3C\/g instead of the less-than character to avoid errors \n        \/\/ in the XML parser.\n        \/\/ Use '\\x26#60;' instead of '<' so that the XML parser\n        \/\/ doesn't go ahead and substitute the less-than character. \n        code_string = code_string.replace(\/\\x3C\/g, '\\x26#60;');\n\n        copyright = 'Copyright 2018 The MathWorks, Inc.';\n\n        w = window.open();\n        d = w.document;\n        d.write('\n\n<pre>\\n');\r\n        d.write(code_string);\r\n\r\n        \/\/ Add copyright line at the bottom if specified.\r\n        if (copyright.length > 0) {\r\n            d.writeln('');\r\n            d.writeln('%%');\r\n            if (copyright.length > 0) {\r\n                d.writeln('% _' + copyright + '_');\r\n            }\r\n        }\r\n\r\n        d.write('<\/pre>\n\n\\n');\n\n        d.title = title + ' (MATLAB code)';\n        d.close();\n    }   \n     --> <\/script><\/p>\n<p style=\"text-align: right; font-size: xx-small; font-weight:lighter;   font-style: italic; color: gray\"><a href=\"javascript:grabCode_c0f6d247caae4ba8a0869999e3f60795()\"><span style=\"font-size: x-small;        font-style: italic;\">Get<br \/>\n      the MATLAB code <noscript>(requires JavaScript)<\/noscript><\/span><\/a><\/p>\n<p>      Published with MATLAB&reg; R2018a<\/p>\n<\/div>\n<p><!--\nc0f6d247caae4ba8a0869999e3f60795 ##### SOURCE BEGIN #####\n%% Path Simplification (and Binary Image Reconstruction!), Made Easy\n%%\n% I've been thinking a lot about how to generate and efficiently store pathsREPLACE_WITH_DASH_DASHand then\n% reconstruct binary images from them. Path simplification is important in\n% many fields, and after some Googling, I found my way to the\n% <https:\/\/en.wikipedia.org\/wiki\/Ramer%E2%80%93Douglas%E2%80%93Peucker_algorithm \"Douglas-Peucker\" algorithm>. (It goes by other names as well, but that\n% seems to be the most common.)\n%\n% It turns out that Douglas-Peucker is used extensively in cartography to\n% simplify the representation of geographical borders. Armed with that\n% information, I came back to MATLAB's doc-search engine and, by typing\n% \"Douglas Peucker,\" I found my way to the <https:\/\/www.mathworks.com\/products\/mapping.html Mapping Toolbox> implementation\n% of |reducem|. From the doc:\n% \n% _The |reducem| function implements a powerful line simplification\n% algorithm (known as Douglas-Peucker) that intelligently selects and deletes visually redundant points._\n%\n% For those of you (us) who have the Mapping Toolbox, you're good to go.\n% Despite the fact that it was implemented for a non-Cartesian geometry, it\n% works nicely, as I will demonstrate. For others, a quick check of the\n% File Exchange reveals _four_ implementations of this particular\n% algorithm:\n%%\n% \n% <<https:\/\/blogs.mathworks.com\/pick\/files\/FexSearch.png>>\n% \n%%\n% Let's create a sample path, and consider each of them (along with |reducem|) in turn:\n%\n\n%% Sample Path\n%\n%  pgon = nsidedpoly(5,'Center',[150 150],'SideLength',100);\n%  verts = pgon.Vertices;\n%  verts = [verts;verts(1,:)];\n%  verts = [verts(1:2:end,:);verts(2:2:end,:)];\n%  mask = poly2mask(verts(:,2),verts(:,1),300,300);\n%  mask = imclose(mask,strel('disk',4));\n%  mask = imfill(mask,'holes');\n%  figure\n%  subplot(1,2,1)\n%  imshow(mask)\n%  title('BW Star')\n%  perim = bwboundaries(mask);\n%  perim = perim{1};\n%  subplot(1,2,2)\n%  plot(perim(:,2),perim(:,1),'b-')\n%  title('Extracted Boundary')\n%  axis equal\n%  set(gca,'xlim',[0,300]+0.5, 'ylim',[0,300]+0.5)\n\n%%\n% \n% <<https:\/\/blogs.mathworks.com\/pick\/files\/star.png>>\n\n%%\n% If you run that code, you'll see that variable _perim_ has 547 [row, col]\n% coordinate pairs. *Most of them are redundant for our purposes*, and we\n% would like to be able to discard them.\n% \n%% Enter the Douglas-Peucker algorithm:\n% \n% * <https:\/\/www.mathworks.com\/matlabcentral\/fileexchange\/61046-douglas-peucker-algorithm Douglas-Peucker Algorithm>, \n% by\n% <https:\/\/www.mathworks.com\/matlabcentral\/profile\/authors\/2664900-reza-ahmadzadeh Reza Ahmadzadeh>:\n% \n%  x = perim(:,2);\n%  y = perim(:,1);\n%  simplified = DouglasPeucker([x';y'],0)\n%\n% <<https:\/\/blogs.mathworks.com\/pick\/files\/DPReza.png>>\n%\n%\n% Wait...what's going on there? Turns out that Reza's implementation considers\n% the distance between the first point and the last. That isn't very useful\n% for a closed curve! I'm going to move on to...\n% \n% * <https:\/\/www.mathworks.com\/matlabcentral\/fileexchange\/36229-douglas-peucker-function Douglas-Peucker Function>, \n% by\n% <https:\/\/www.mathworks.com\/matlabcentral\/profile\/authors\/3289512-ranadeer Ranadeer>:\n% \n%  [list_pts,pnts]= douglas_peucker([x,y], 0);\n%  size(pnts)\n%  ans =\n%       2     2\n%\n% Here again, the results match Reza's (the two returned points are the\n% point on the \"upper left\" prong, and are coincident). Ranadeer's\n% implementation also relies on distances from the endpoints of the curve. (Note\n% that changing the tolerance in either Reza's or Ranadeer's implementation\n% has no effect.) Let's try...\n%\n% \n% * <https:\/\/www.mathworks.com\/matlabcentral\/fileexchange\/41986-ramer-douglas-peucker-algorithm-demo Ramer-Douglas-Peucker>, \n% by\n% <https:\/\/www.mathworks.com\/matlabcentral\/profile\/authors\/4185006-ligong-han Ligong Han>:\n% \n%  simplified = DouglasPeucker([x,y],10);\n%\n% <<https:\/\/blogs.mathworks.com\/pick\/files\/DPHan.png>>\n%\n% Beautiful! That's what I'm talking about! Conveniently, Ligong coded to\n% allow an additional input flag to show the data:\n%\n% <<https:\/\/blogs.mathworks.com\/pick\/files\/DPHan2.png>>\n%\n% Nice!\n%\n% Next up, <https:\/\/www.mathworks.com\/matlabcentral\/fileexchange\/21132-line-simplification Line Simplification>,\n% by\n% <https:\/\/www.mathworks.com\/matlabcentral\/profile\/authors\/870595-wolfgang-schwanghart Wolfgang Schwanghart>:\n%\n%  [ps,ix] = dpsimplify([x,y],10);\n%\n% <<https:\/\/blogs.mathworks.com\/pick\/files\/DPSchwanghart.png>>\n%\n% Again, that's lovely! (Note that Wolfgang also outputs idxREPLACE_WITH_DASH_DASHan index of\n% the points that are retained. That could be useful!)\n% \n% Finally, consider the Mapping Toolbox implementation. Even though it is\n% implemented for geospatial (latitude\/longitude) data, it works\n% beautifully!\n%\n% * |reducem|:\n%\n%  [xReduced,yReduced] = reducem(perim(:,2),perim(:,1),10);\n%  plot(xReduced,yReduced,'r.','Markersize',14)\n%\n% <<https:\/\/blogs.mathworks.com\/pick\/files\/reducem.png>>\n%\n\n%% Mask reconstruction\n% So clearly we have some great options, even for simplifying closed\n% curves. Once you've got that reduced set of vertices, note that you can\n% readily reconstruct the polygon using poly2mask():\n%\n% mask2 = poly2mask(ps(:,1),ps(:,2), 300, 300);\n%\n% <<https:\/\/blogs.mathworks.com\/pick\/files\/Reconstructed.png>>\n%\n% It is important to note that, while those binary images are nearly\n% identical, they are not _exactly_ identical. Why? Steve (in one of his\n% stellar \n% <https:\/\/blogs.mathworks.com\/steve\/2014\/03\/27\/comparing-the-geometries-of-bwboundaries-and-poly2mask\/ blog posts> on\n% image processing) discusses that issue. It turns out that we _can_ get\n% identical results, if we upsample the path, with nearest-neighbor\n% interpolation, by a factor of 3!\n%%\n% As always, I welcome your\n% <http:\/\/blogs.mathworks.com\/pick\/?p=10035#respond thoughts and comments>.\n##### SOURCE END ##### c0f6d247caae4ba8a0869999e3f60795\n--><\/p>\n","protected":false},"excerpt":{"rendered":"<div class=\"overview-image\"><img src=\"https:\/\/blogs.mathworks.com\/pick\/files\/star.png\" class=\"img-responsive attachment-post-thumbnail size-post-thumbnail wp-post-image\" alt=\"\" decoding=\"async\" loading=\"lazy\" \/><\/div>\n<p>\nContents<\/p>\n<p>Sample Path<br \/>\nEnter the Douglas-Peucker algorithm:<br \/>\nMask reconstruction<\/p>\n<p>I&#8217;ve been thinking a lot about how to generate and efficiently store paths&#8211;and then reconstruct binary&#8230; <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/pick\/2018\/09\/07\/path-simplification-and-binary-image-reconstruction-made-easy\/\">read more >><\/a><\/p>\n","protected":false},"author":34,"featured_media":10051,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[16],"tags":[],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/posts\/10035"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/users\/34"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/comments?post=10035"}],"version-history":[{"count":8,"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/posts\/10035\/revisions"}],"predecessor-version":[{"id":10093,"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/posts\/10035\/revisions\/10093"}],"wp:featuredmedia":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/media\/10051"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/media?parent=10035"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/categories?post=10035"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/pick\/wp-json\/wp\/v2\/tags?post=10035"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}