{"id":809,"date":"2013-03-21T15:14:57","date_gmt":"2013-03-21T19:14:57","guid":{"rendered":"https:\/\/blogs.mathworks.com\/steve\/?p=809"},"modified":"2019-11-01T09:13:27","modified_gmt":"2019-11-01T13:13:27","slug":"revisiting-dctdemo-part-3","status":"publish","type":"post","link":"https:\/\/blogs.mathworks.com\/steve\/2013\/03\/21\/revisiting-dctdemo-part-3\/","title":{"rendered":"Revisiting dctdemo &#8211; part 3"},"content":{"rendered":"<div class=\"content\"><p>This is the third part of my plan to rewrite an Image Processing Toolbox example from 20 years ago using more modern MATLAB language features. I got the idea from Dave Garrison's <a href=\"https:\/\/www.mathworks.com\/company\/newsletters\/articles\/writing-apps-in-matlab.html\">recent article<\/a> on writing MATLAB apps.<\/p><p>Here's the old app I'm trying to reinvent:<\/p><p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dctdemo-screen-shot.png\" alt=\"\"> <\/p><p>And here's what I had working the last time:<\/p><p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dct-compression-v2.png\" alt=\"\"> <\/p><p>The code isn't doing any actual DCT computations yet. The idea is to compute 8-by-8 block DCTs, zero out the least significant DCT coefficients of each block, and then reconstruct the image by computing the inverse DCT of each 8-by-8 block.<\/p><p>Let me show how this will work on the sample image I'm using for the app:<\/p><pre class=\"codeinput\">pout = imread(<span class=\"string\">'pout.tif'<\/span>);\r\npout2 = pout(1:240,:);\r\nI = im2double(adapthisteq(pout2));\r\nimshow(I)\r\n<\/pre><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dctdemo_part3_01.png\" alt=\"\"> <p>First, use <tt>blockproc<\/tt> to compute the 8-by-8 2-D DCTs.<\/p><pre class=\"codeinput\">f = @(block) dct2(block.data);\r\nA = blockproc(I,[8 8],f);\r\nimshow(A,[])\r\ntitle(<span class=\"string\">'DCT coefficients'<\/span>)\r\n<\/pre><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dctdemo_part3_02.png\" alt=\"\"> <p>Second, I want to compute the variances of each of the 64 coefficients. That is, compute the variance of the (1,1) DCT coefficients from each block, the variance of the (1,2) DCT coefficients from each block, and so on. I'll use the computed variances to determine which DCT coefficients to keep and which to zero out. To do this, I'll rearrange the elements from the previous step so that the corresponding DCT coefficients are together in the columns of a 64-column matrix.<\/p><pre class=\"codeinput\">B = im2col(A,[8 8],<span class=\"string\">'distinct'<\/span>)';\r\nvars = var(B);\r\nplot(vars)\r\ntitle(<span class=\"string\">'Variances of the 64 DCT coefficients'<\/span>)\r\n<\/pre><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dctdemo_part3_03.png\" alt=\"\"> <p>Next, suppose we want to keep the 3 block DCT coefficients with the highest variances. I'll use <tt>sort<\/tt> to figure out which ones to keep.<\/p><pre class=\"codeinput\">[~,idx] = sort(vars,<span class=\"string\">'descend'<\/span>);\r\nkeep = idx(1:3)\r\n<\/pre><pre class=\"codeoutput\">\r\nkeep =\r\n\r\n     1     9     2\r\n\r\n<\/pre><p>For an 8-by-8 block of DCT coefficients, index values 1, 9, and 2 correspond to the upper-left DCT coefficient and the two coefficients just the right and down from it. Let's keep just those coefficients.<\/p><pre class=\"codeinput\">B2 = zeros(size(B));\r\nB2(:,keep) = B(:,keep);\r\n<\/pre><p>Now I'll rearrange the columns back into blocks and reconstruct the image by computing the inverse DCT of each block.<\/p><pre class=\"codeinput\">C = col2im(B2',[8 8],size(I),<span class=\"string\">'distinct'<\/span>);\r\nfinv = @(block) idct2(block.data);\r\nD = blockproc(C,[8 8],finv);\r\nimshow(D)\r\ntitle(<span class=\"string\">'Image reconstructed from 3 DCT coefficients per block'<\/span>)\r\n<\/pre><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dctdemo_part3_04.png\" alt=\"\"> <p>This algorithm code into a local function to go at the bottom of my app code file, like this:<\/p><pre>function I2 = reconstructImage(I,n)\r\n% Reconstruct the image from n of the DCT coefficients in each 8-by-8\r\n% block. Select the n coefficients with the largest variance across the\r\n% image.<\/pre><pre>% Compute 8-by-8 block DCTs.\r\nf = @(block) dct2(block.data);\r\nA = blockproc(I,[8 8],f);<\/pre><pre>% Compute DCT coefficient variances and decide\r\n% which to keep.\r\nB = im2col(A,[8 8],'distinct')';\r\nvars = var(B);\r\n[~,idx] = sort(vars,'descend');\r\nkeep = idx(1:n);<\/pre><pre>% Zero out the DCT coefficients we are not keeping.\r\nB2 = zeros(size(B));\r\nB2(:,keep) = B(:,keep);<\/pre><pre>% Reconstruct image using 8-by-8 block inverse\r\n% DCTs.\r\nC = col2im(B2',[8 8],size(I),'distinct');\r\nfinv = @(block) idct2(block.data);\r\nI2 = blockproc(C,[8 8],finv);\r\nend<\/pre><p>I already have an <tt>update<\/tt> method in my class, so I can add code to that method that will compute the reconstructed image and error image. Here's the revised <tt>update<\/tt> method:<\/p><pre>function update(app)\r\n    recon_image = reconstructImage(app.OriginalImage, ...\r\n        app.NumDCTCoefficients);<\/pre><pre>    diff_image = imabsdiff(app.OriginalImage, recon_image);<\/pre><pre>    imshow(app.OriginalImage,'Parent',app.OriginalImageAxes);\r\n    title(app.OriginalImageAxes,'Original Image');<\/pre><pre>    imshow(recon_image,'Parent',app.ReconstructedImageAxes);\r\n    title(app.ReconstructedImageAxes,'Reconstructed Image');<\/pre><pre>    imshow(diff_image,[],'Parent',app.ErrorImageAxes);\r\n    title(app.ErrorImageAxes,'Error Image');<\/pre><pre>    imshow(zeros(size(app.OriginalImage)),'Parent',app.MaskAxes);\r\n    title(app.MaskAxes,'DCT Coefficient Mask');\r\n    drawnow;\r\nend<\/pre><p>Although I don't have all the interactive behavior wired up yet, and don't have the DCT coefficient mask visualization done yet, the app runs and looks a bit more finished than it did before. It now shows the results for reconstructing the image from 3 DCT coefficients per block.<\/p><p><img decoding=\"async\" vspace=\"5\" hspace=\"5\" src=\"https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dct-compression-example-v3.png\" alt=\"\"> <\/p><p>I think it'll take just one more post to finish making the rest of this app work.<\/p>\r\n<div><ul><li><a href=\"https:\/\/blogs.mathworks.com\/steve\/2013\/02\/08\/revisiting-dctdemo-part-1\/\">Revisiting dctdemo - part 1<\/a><\/li><li><a href=\"https:\/\/blogs.mathworks.com\/steve\/2013\/02\/21\/revisiting-dctdemo-part-2\/\">Revisiting dctdemo - part 2<\/a><\/li><li><a href=\"https:\/\/blogs.mathworks.com\/steve\/2013\/03\/21\/revisiting-dctdemo-part-3\/\">Revisiting dctdemo - part 3<\/a><\/li><li><a href=\"https:\/\/blogs.mathworks.com\/steve\/2013\/04\/12\/revisiting-dctdemo-part-4\/\">Revisiting dctdemo - part 4<\/a><\/li><\/ul><\/div>\r\n<script language=\"JavaScript\"> <!-- \r\n    function grabCode_03724e04a81d4ac29734fd49392c5159() {\r\n        \/\/ Remember the title so we can use it in the new page\r\n        title = document.title;\r\n\r\n        \/\/ Break up these strings so that their presence\r\n        \/\/ in the Javascript doesn't mess up the search for\r\n        \/\/ the MATLAB code.\r\n        t1='03724e04a81d4ac29734fd49392c5159 ' + '##### ' + 'SOURCE BEGIN' + ' #####';\r\n        t2='##### ' + 'SOURCE END' + ' #####' + ' 03724e04a81d4ac29734fd49392c5159';\r\n    \r\n        b=document.getElementsByTagName('body')[0];\r\n        i1=b.innerHTML.indexOf(t1)+t1.length;\r\n        i2=b.innerHTML.indexOf(t2);\r\n \r\n        code_string = b.innerHTML.substring(i1, i2);\r\n        code_string = code_string.replace(\/REPLACE_WITH_DASH_DASH\/g,'--');\r\n\r\n        \/\/ Use \/x3C\/g instead of the less-than character to avoid errors \r\n        \/\/ in the XML parser.\r\n        \/\/ Use '\\x26#60;' instead of '<' so that the XML parser\r\n        \/\/ doesn't go ahead and substitute the less-than character. \r\n        code_string = code_string.replace(\/\\x3C\/g, '\\x26#60;');\r\n\r\n        copyright = 'Copyright 2013 The MathWorks, Inc.';\r\n\r\n        w = window.open();\r\n        d = w.document;\r\n        d.write('<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');\r\n\r\n        d.title = title + ' (MATLAB code)';\r\n        d.close();\r\n    }   \r\n     --> <\/script><p style=\"text-align: right; font-size: xx-small; font-weight:lighter;   font-style: italic; color: gray\"><br><a href=\"javascript:grabCode_03724e04a81d4ac29734fd49392c5159()\"><span style=\"font-size: x-small;        font-style: italic;\">Get \r\n      the MATLAB code <noscript>(requires JavaScript)<\/noscript><\/span><\/a><br><br>\r\n      Published with MATLAB&reg; R2013a<br><\/p><p class=\"footer\"><br>\r\n      Published with MATLAB&reg; R2013a<br><\/p><\/div><!--\r\n03724e04a81d4ac29734fd49392c5159 ##### SOURCE BEGIN #####\r\n%%\r\n% This is the third part of my plan to rewrite an Image Processing Toolbox\r\n% example from 20 years ago using more modern MATLAB language features. I\r\n% got the idea from Dave Garrison's \r\n% <https:\/\/www.mathworks.com\/company\/newsletters\/articles\/writing-apps-in-matlab.html \r\n% recent article> on writing MATLAB apps.\r\n%\r\n% Here's the old app I'm trying to reinvent:\r\n%\r\n% <<https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dctdemo-screen-shot.png>>\r\n%\r\n% And here's what I had working the last time:\r\n%\r\n% <<https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dct-compression-v2.png>>\r\n%\r\n% The code isn't doing any actual DCT computations yet. The idea is to\r\n% compute 8-by-8 block DCTs, zero out the least significant DCT\r\n% coefficients of each block, and then reconstruct the image by computing\r\n% the inverse DCT of each 8-by-8 block.\r\n%\r\n% Let me show how this will work on the sample image I'm using for the app:\r\n\r\npout = imread('pout.tif');\r\npout2 = pout(1:240,:);\r\nI = im2double(adapthisteq(pout2));\r\nimshow(I)\r\n\r\n%%\r\n% First, use |blockproc| to compute the 8-by-8 2-D DCTs.\r\n\r\nf = @(block) dct2(block.data);\r\nA = blockproc(I,[8 8],f);\r\nimshow(A,[])\r\ntitle('DCT coefficients')\r\n\r\n%%\r\n% Second, I want to compute the variances of each of the 64 coefficients.\r\n% That is, compute the variance of the (1,1) DCT coefficients from each\r\n% block, the variance of the (1,2) DCT coefficients from each block, and so\r\n% on. I'll use the computed variances to determine which DCT coefficients\r\n% to keep and which to zero out. To do this, I'll rearrange the elements\r\n% from the previous step so that the corresponding DCT coefficients are\r\n% together in the columns of a 64-column matrix.\r\n\r\nB = im2col(A,[8 8],'distinct')';\r\nvars = var(B);\r\nplot(vars)\r\ntitle('Variances of the 64 DCT coefficients')\r\n\r\n%%\r\n% Next, suppose we want to keep the 3 block DCT coefficients with the\r\n% highest variances. I'll use |sort| to figure out which ones to keep.\r\n\r\n[~,idx] = sort(vars,'descend');\r\nkeep = idx(1:3)\r\n\r\n%%\r\n% For an 8-by-8 block of DCT coefficients, index values 1, 9, and 2\r\n% correspond to the upper-left DCT coefficient and the two coefficients\r\n% just the right and down from it. Let's keep just those coefficients.\r\n\r\nB2 = zeros(size(B));\r\nB2(:,keep) = B(:,keep);\r\n\r\n%%\r\n% Now I'll rearrange the columns back into blocks and reconstruct the image\r\n% by computing the inverse DCT of each block.\r\n\r\nC = col2im(B2',[8 8],size(I),'distinct');\r\nfinv = @(block) idct2(block.data);\r\nD = blockproc(C,[8 8],finv);\r\nimshow(D)\r\ntitle('Image reconstructed from 3 DCT coefficients per block')\r\n\r\n%%\r\n% This algorithm code into a local function to go at the\r\n% bottom of my app code file, like this:\r\n%\r\n%  function I2 = reconstructImage(I,n)\r\n%  % Reconstruct the image from n of the DCT coefficients in each 8-by-8\r\n%  % block. Select the n coefficients with the largest variance across the\r\n%  % image.\r\n%  \r\n%  % Compute 8-by-8 block DCTs.\r\n%  f = @(block) dct2(block.data);\r\n%  A = blockproc(I,[8 8],f);\r\n%  \r\n%  % Compute DCT coefficient variances and decide\r\n%  % which to keep.\r\n%  B = im2col(A,[8 8],'distinct')';\r\n%  vars = var(B);\r\n%  [~,idx] = sort(vars,'descend');\r\n%  keep = idx(1:n);\r\n%  \r\n%  % Zero out the DCT coefficients we are not keeping.\r\n%  B2 = zeros(size(B));\r\n%  B2(:,keep) = B(:,keep);\r\n%  \r\n%  % Reconstruct image using 8-by-8 block inverse\r\n%  % DCTs.\r\n%  C = col2im(B2',[8 8],size(I),'distinct');\r\n%  finv = @(block) idct2(block.data);\r\n%  I2 = blockproc(C,[8 8],finv);\r\n%  end\r\n%\r\n% I already have an |update| method in my class, so I can add code to that\r\n% method that will compute the reconstructed image and error image. Here's\r\n% the revised |update| method:\r\n%\r\n%  function update(app)\r\n%      recon_image = reconstructImage(app.OriginalImage, ...\r\n%          app.NumDCTCoefficients);\r\n%      \r\n%      diff_image = imabsdiff(app.OriginalImage, recon_image);\r\n%      \r\n%      imshow(app.OriginalImage,'Parent',app.OriginalImageAxes);\r\n%      title(app.OriginalImageAxes,'Original Image');\r\n%      \r\n%      imshow(recon_image,'Parent',app.ReconstructedImageAxes);\r\n%      title(app.ReconstructedImageAxes,'Reconstructed Image');\r\n%      \r\n%      imshow(diff_image,[],'Parent',app.ErrorImageAxes);\r\n%      title(app.ErrorImageAxes,'Error Image');\r\n%      \r\n%      imshow(zeros(size(app.OriginalImage)),'Parent',app.MaskAxes);\r\n%      title(app.MaskAxes,'DCT Coefficient Mask');\r\n%      drawnow;\r\n%  end\r\n%\r\n% Although I don't have all the interactive behavior wired up yet, and\r\n% don't have the DCT coefficient mask visualization done yet, the app runs\r\n% and looks a bit more finished than it did before. It now shows the\r\n% results for reconstructing the image from 3 DCT coefficients per block.\r\n%\r\n% <<https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dct-compression-example-v3.png>>\r\n%\r\n% I think it'll take just one more post to finish making the rest of this\r\n% app work.\r\n\r\n##### SOURCE END ##### 03724e04a81d4ac29734fd49392c5159\r\n-->","protected":false},"excerpt":{"rendered":"<div class=\"overview-image\"><img decoding=\"async\"  class=\"img-responsive\" src=\"https:\/\/blogs.mathworks.com\/images\/steve\/2013\/dct-compression-example-v3.png\" onError=\"this.style.display ='none';\" \/><\/div><p>This is the third part of my plan to rewrite an Image Processing Toolbox example from 20 years ago using more modern MATLAB language features. I got the idea from Dave Garrison's recent article on... <a class=\"read-more\" href=\"https:\/\/blogs.mathworks.com\/steve\/2013\/03\/21\/revisiting-dctdemo-part-3\/\">read more >><\/a><\/p>","protected":false},"author":42,"featured_media":0,"comment_status":"open","ping_status":"closed","sticky":false,"template":"","format":"standard","meta":[],"categories":[1],"tags":[134,711,989,204,981,991,987,390,993,76,36,68,190,484,52,432,130],"_links":{"self":[{"href":"https:\/\/blogs.mathworks.com\/steve\/wp-json\/wp\/v2\/posts\/809"}],"collection":[{"href":"https:\/\/blogs.mathworks.com\/steve\/wp-json\/wp\/v2\/posts"}],"about":[{"href":"https:\/\/blogs.mathworks.com\/steve\/wp-json\/wp\/v2\/types\/post"}],"author":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/steve\/wp-json\/wp\/v2\/users\/42"}],"replies":[{"embeddable":true,"href":"https:\/\/blogs.mathworks.com\/steve\/wp-json\/wp\/v2\/comments?post=809"}],"version-history":[{"count":5,"href":"https:\/\/blogs.mathworks.com\/steve\/wp-json\/wp\/v2\/posts\/809\/revisions"}],"predecessor-version":[{"id":3815,"href":"https:\/\/blogs.mathworks.com\/steve\/wp-json\/wp\/v2\/posts\/809\/revisions\/3815"}],"wp:attachment":[{"href":"https:\/\/blogs.mathworks.com\/steve\/wp-json\/wp\/v2\/media?parent=809"}],"wp:term":[{"taxonomy":"category","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/steve\/wp-json\/wp\/v2\/categories?post=809"},{"taxonomy":"post_tag","embeddable":true,"href":"https:\/\/blogs.mathworks.com\/steve\/wp-json\/wp\/v2\/tags?post=809"}],"curies":[{"name":"wp","href":"https:\/\/api.w.org\/{rel}","templated":true}]}}