| idamnjanovic@6 | 1 %% DICTIONARY LEARNING FOR IMAGE DENOISING | 
| ivan@107 | 2 | 
| idamnjanovic@25 | 3 % | 
| idamnjanovic@6 | 4 %   This file contains an example of how SMALLbox can be used to test different | 
| idamnjanovic@6 | 5 %   dictionary learning techniques in Image Denoising problem. | 
| idamnjanovic@6 | 6 %   It calls generateImageDenoiseProblem that will let you to choose image, | 
| idamnjanovic@6 | 7 %   add noise and use noisy image to generate training set for dictionary | 
| idamnjanovic@6 | 8 %   learning. | 
| idamnjanovic@6 | 9 %   Three dictionary learning techniques were compared: | 
| idamnjanovic@6 | 10 %   -   KSVD - M. Elad, R. Rubinstein, and M. Zibulevsky, "Efficient | 
| idamnjanovic@6 | 11 %              Implementation of the K-SVD Algorithm using Batch Orthogonal | 
| idamnjanovic@6 | 12 %              Matching Pursuit", Technical Report - CS, Technion, April 2008. | 
| idamnjanovic@6 | 13 %   -   KSVDS - R. Rubinstein, M. Zibulevsky, and M. Elad, "Learning Sparse | 
| idamnjanovic@6 | 14 %               Dictionaries for Sparse Signal Approximation", Technical | 
| idamnjanovic@6 | 15 %               Report - CS, Technion, June 2009. | 
| idamnjanovic@6 | 16 %   -   SPAMS - J. Mairal, F. Bach, J. Ponce and G. Sapiro. Online | 
| idamnjanovic@6 | 17 %               Dictionary Learning for Sparse Coding. International | 
| idamnjanovic@6 | 18 %               Conference on Machine Learning,Montreal, Canada, 2009 | 
| idamnjanovic@6 | 19 % | 
| ivan@107 | 20 | 
| ivan@107 | 21 % | 
| ivan@107 | 22 %   Centre for Digital Music, Queen Mary, University of London. | 
| ivan@107 | 23 %   This file copyright 2009 Ivan Damnjanovic. | 
| ivan@107 | 24 % | 
| ivan@107 | 25 %   This program is free software; you can redistribute it and/or | 
| ivan@107 | 26 %   modify it under the terms of the GNU General Public License as | 
| ivan@107 | 27 %   published by the Free Software Foundation; either version 2 of the | 
| ivan@107 | 28 %   License, or (at your option) any later version.  See the file | 
| ivan@107 | 29 %   COPYING included with this distribution for more information. | 
| idamnjanovic@6 | 30 % | 
| idamnjanovic@6 | 31 %% | 
| idamnjanovic@6 | 32 | 
| idamnjanovic@6 | 33 clear; | 
| idamnjanovic@6 | 34 | 
| idamnjanovic@6 | 35 %   If you want to load the image outside of generateImageDenoiseProblem | 
| idamnjanovic@6 | 36 %   function uncomment following lines. This can be useful if you want to | 
| idamnjanovic@6 | 37 %   denoise more then one image for example. | 
| idamnjanovic@6 | 38 | 
| idamnjanovic@6 | 39 % TMPpath=pwd; | 
| idamnjanovic@6 | 40 % FS=filesep; | 
| idamnjanovic@6 | 41 % [pathstr1, name, ext, versn] = fileparts(which('SMALLboxSetup.m')); | 
| idamnjanovic@6 | 42 % cd([pathstr1,FS,'data',FS,'images']); | 
| idamnjanovic@6 | 43 % [filename,pathname] = uigetfile({'*.png;'},'Select a file containin pre-calculated notes'); | 
| idamnjanovic@6 | 44 % [pathstr, name, ext, versn] = fileparts(filename); | 
| idamnjanovic@6 | 45 % test_image = imread(filename); | 
| idamnjanovic@6 | 46 % test_image = double(test_image); | 
| idamnjanovic@6 | 47 % cd(TMPpath); | 
| idamnjanovic@6 | 48 % SMALL.Problem.name=name; | 
| idamnjanovic@6 | 49 | 
| idamnjanovic@6 | 50 | 
| idamnjanovic@6 | 51 % Defining Image Denoising Problem as Dictionary Learning | 
| idamnjanovic@6 | 52 % Problem. As an input we set the number of training patches. | 
| idamnjanovic@6 | 53 | 
| idamnjanovic@6 | 54 SMALL.Problem = generateImageDenoiseProblem('', 40000); | 
| idamnjanovic@6 | 55 | 
| idamnjanovic@6 | 56 | 
| idamnjanovic@6 | 57 %% | 
| idamnjanovic@6 | 58 %   Use KSVD Dictionary Learning Algorithm to Learn overcomplete dictionary | 
| idamnjanovic@6 | 59 | 
| idamnjanovic@6 | 60 %   Initialising Dictionary structure | 
| idamnjanovic@6 | 61 %   Setting Dictionary structure fields (toolbox, name, param, D and time) | 
| idamnjanovic@6 | 62 %   to zero values | 
| idamnjanovic@6 | 63 | 
| idamnjanovic@6 | 64 SMALL.DL(1)=SMALL_init_DL(); | 
| idamnjanovic@6 | 65 | 
| idamnjanovic@6 | 66 % Defining the parameters needed for dictionary learning | 
| idamnjanovic@6 | 67 | 
| idamnjanovic@6 | 68 SMALL.DL(1).toolbox = 'KSVD'; | 
| idamnjanovic@6 | 69 SMALL.DL(1).name = 'ksvd'; | 
| idamnjanovic@6 | 70 | 
| idamnjanovic@6 | 71 %   Defining the parameters for KSVD | 
| idamnjanovic@6 | 72 %   In this example we are learning 256 atoms in 20 iterations, so that | 
| idamnjanovic@6 | 73 %   every patch in the training set can be represented with target error in | 
| idamnjanovic@6 | 74 %   L2-norm (EData) | 
| idamnjanovic@6 | 75 %   Type help ksvd in MATLAB prompt for more options. | 
| idamnjanovic@6 | 76 | 
| idamnjanovic@6 | 77 Edata=sqrt(prod(SMALL.Problem.blocksize)) * SMALL.Problem.sigma * SMALL.Problem.gain; | 
| ivan@107 | 78 maxatoms = floor(prod(SMALL.Problem.blocksize)/2); | 
| ivan@107 | 79 | 
| idamnjanovic@6 | 80 SMALL.DL(1).param=struct(... | 
| idamnjanovic@6 | 81     'Edata', Edata,... | 
| idamnjanovic@6 | 82     'initdict', SMALL.Problem.initdict,... | 
| idamnjanovic@6 | 83     'dictsize', SMALL.Problem.p,... | 
| idamnjanovic@6 | 84     'iternum', 20,... | 
| idamnjanovic@6 | 85     'memusage', 'high'); | 
| idamnjanovic@6 | 86 | 
| idamnjanovic@6 | 87 %   Learn the dictionary | 
| idamnjanovic@6 | 88 | 
| idamnjanovic@6 | 89 SMALL.DL(1) = SMALL_learn(SMALL.Problem, SMALL.DL(1)); | 
| idamnjanovic@6 | 90 | 
| idamnjanovic@6 | 91 %   Set SMALL.Problem.A dictionary | 
| idamnjanovic@6 | 92 %   (backward compatiblity with SPARCO: solver structure communicate | 
| idamnjanovic@6 | 93 %   only with Problem structure, ie no direct communication between DL and | 
| idamnjanovic@6 | 94 %   solver structures) | 
| idamnjanovic@6 | 95 | 
| idamnjanovic@6 | 96 SMALL.Problem.A = SMALL.DL(1).D; | 
| ivan@107 | 97 SMALL.Problem.reconstruct = @(x) ImgDenoise_reconstruct(x, SMALL.Problem); | 
| idamnjanovic@6 | 98 | 
| idamnjanovic@6 | 99 %% | 
| idamnjanovic@6 | 100 %   Initialising solver structure | 
| idamnjanovic@6 | 101 %   Setting solver structure fields (toolbox, name, param, solution, | 
| idamnjanovic@6 | 102 %   reconstructed and time) to zero values | 
| idamnjanovic@6 | 103 | 
| idamnjanovic@6 | 104 SMALL.solver(1)=SMALL_init_solver; | 
| idamnjanovic@6 | 105 | 
| idamnjanovic@6 | 106 % Defining the parameters needed for image denoising | 
| idamnjanovic@6 | 107 | 
| idamnjanovic@6 | 108 SMALL.solver(1).toolbox='ompbox'; | 
| ivan@107 | 109 SMALL.solver(1).name='omp2'; | 
| ivan@107 | 110 SMALL.solver(1).param=struct(... | 
| ivan@107 | 111     'epsilon',Edata,... | 
| ivan@107 | 112     'maxatoms', maxatoms); | 
| idamnjanovic@6 | 113 | 
| ivan@107 | 114 %   Denoising the image - find the sparse solution in the learned | 
| ivan@107 | 115 %   dictionary for all patches in the image and the end it uses | 
| ivan@107 | 116 %   reconstruction function to reconstruct the patches and put them into a | 
| ivan@107 | 117 %   denoised image | 
| idamnjanovic@6 | 118 | 
| ivan@107 | 119 SMALL.solver(1)=SMALL_solve(SMALL.Problem, SMALL.solver(1)); | 
| ivan@107 | 120 | 
| ivan@107 | 121 %   Show PSNR after reconstruction | 
| ivan@107 | 122 | 
| ivan@107 | 123 SMALL.solver(1).reconstructed.psnr | 
| idamnjanovic@6 | 124 | 
| idamnjanovic@6 | 125 %% | 
| idamnjanovic@6 | 126 % Use KSVDS Dictionary Learning Algorithm to denoise image | 
| idamnjanovic@6 | 127 | 
| idamnjanovic@6 | 128 %   Initialising solver structure | 
| idamnjanovic@6 | 129 %   Setting solver structure fields (toolbox, name, param, solution, | 
| idamnjanovic@6 | 130 %   reconstructed and time) to zero values | 
| idamnjanovic@6 | 131 | 
| idamnjanovic@6 | 132 SMALL.DL(2)=SMALL_init_DL(); | 
| idamnjanovic@6 | 133 | 
| idamnjanovic@6 | 134 % Defining the parameters needed for dictionary learning | 
| idamnjanovic@6 | 135 | 
| idamnjanovic@6 | 136 SMALL.DL(2).toolbox = 'KSVDS'; | 
| idamnjanovic@6 | 137 SMALL.DL(2).name = 'ksvds'; | 
| idamnjanovic@6 | 138 | 
| idamnjanovic@6 | 139 %   Defining the parameters for KSVDS | 
| idamnjanovic@6 | 140 %   In this example we are learning 256 atoms in 20 iterations, so that | 
| idamnjanovic@6 | 141 %   every patch in the training set can be represented with target error in | 
| idamnjanovic@6 | 142 %   L2-norm (EDataS). We also impose "double sparsity" - dictionary itself | 
| idamnjanovic@6 | 143 %   has to be sparse in the given base dictionary (Tdict - number of | 
| idamnjanovic@6 | 144 %   nonzero elements per atom). | 
| idamnjanovic@6 | 145 %   Type help ksvds in MATLAB prompt for more options. | 
| idamnjanovic@6 | 146 | 
| idamnjanovic@6 | 147 EdataS=sqrt(prod(SMALL.Problem.blocksize)) * SMALL.Problem.sigma * SMALL.Problem.gain; | 
| idamnjanovic@6 | 148 SMALL.DL(2).param=struct(... | 
| idamnjanovic@6 | 149     'Edata', EdataS, ... | 
| idamnjanovic@6 | 150     'Tdict', 6,... | 
| idamnjanovic@6 | 151     'stepsize', 1,... | 
| idamnjanovic@6 | 152     'dictsize', SMALL.Problem.p,... | 
| idamnjanovic@6 | 153     'iternum', 20,... | 
| idamnjanovic@6 | 154     'memusage', 'high'); | 
| idamnjanovic@6 | 155 SMALL.DL(2).param.initA = speye(SMALL.Problem.p); | 
| idamnjanovic@6 | 156 SMALL.DL(2).param.basedict{1} = odctdict(8,16); | 
| idamnjanovic@6 | 157 SMALL.DL(2).param.basedict{2} = odctdict(8,16); | 
| idamnjanovic@6 | 158 | 
| idamnjanovic@6 | 159 % Learn the dictionary | 
| idamnjanovic@6 | 160 | 
| idamnjanovic@6 | 161 SMALL.DL(2) = SMALL_learn(SMALL.Problem, SMALL.DL(2)); | 
| idamnjanovic@6 | 162 | 
| idamnjanovic@6 | 163 %   Set SMALL.Problem.A dictionary and SMALL.Problem.basedictionary | 
| idamnjanovic@6 | 164 %   (backward compatiblity with SPARCO: solver structure communicate | 
| idamnjanovic@6 | 165 %   only with Problem structure, ie no direct communication between DL and | 
| idamnjanovic@6 | 166 %   solver structures) | 
| idamnjanovic@6 | 167 | 
| idamnjanovic@6 | 168 SMALL.Problem.A = SMALL.DL(2).D; | 
| idamnjanovic@6 | 169 SMALL.Problem.basedict{1} = SMALL.DL(2).param.basedict{1}; | 
| idamnjanovic@6 | 170 SMALL.Problem.basedict{2} = SMALL.DL(2).param.basedict{2}; | 
| idamnjanovic@6 | 171 | 
| ivan@107 | 172 %   Setting up reconstruction function | 
| ivan@107 | 173 | 
| ivan@107 | 174 SparseDict=1; | 
| ivan@107 | 175 SMALL.Problem.reconstruct = @(x) ImgDenoise_reconstruct(x, SMALL.Problem, SparseDict); | 
| ivan@107 | 176 | 
| idamnjanovic@6 | 177 %   Initialising solver structure | 
| idamnjanovic@6 | 178 %   Setting solver structure fields (toolbox, name, param, solution, | 
| idamnjanovic@6 | 179 %   reconstructed and time) to zero values | 
| idamnjanovic@6 | 180 | 
| idamnjanovic@6 | 181 SMALL.solver(2)=SMALL_init_solver; | 
| idamnjanovic@6 | 182 | 
| idamnjanovic@6 | 183 % Defining the parameters needed for image denoising | 
| idamnjanovic@6 | 184 | 
| idamnjanovic@6 | 185 SMALL.solver(2).toolbox='ompsbox'; | 
| ivan@107 | 186 SMALL.solver(2).name='omps2'; | 
| ivan@107 | 187 SMALL.solver(2).param=struct(... | 
| ivan@107 | 188     'epsilon',Edata,... | 
| ivan@107 | 189     'maxatoms', maxatoms); | 
| idamnjanovic@6 | 190 | 
| ivan@107 | 191 %   Denoising the image - find the sparse solution in the learned | 
| ivan@107 | 192 %   dictionary for all patches in the image and the end it uses | 
| ivan@107 | 193 %   reconstruction function to reconstruct the patches and put them into a | 
| ivan@107 | 194 %   denoised image | 
| idamnjanovic@6 | 195 | 
| ivan@107 | 196 SMALL.solver(2)=SMALL_solve(SMALL.Problem, SMALL.solver(2)); | 
| idamnjanovic@6 | 197 | 
| ivan@107 | 198 %% | 
| ivan@107 | 199 %   Use SPAMS Online Dictionary Learning Algorithm | 
| ivan@107 | 200 %   to Learn overcomplete dictionary (Julien Mairal 2009) | 
| ivan@107 | 201 %   (If you have not installed SPAMS please comment the following two cells) | 
| ivan@107 | 202 | 
| ivan@107 | 203 %   Initialising Dictionary structure | 
| ivan@107 | 204 %   Setting Dictionary structure fields (toolbox, name, param, D and time) | 
| ivan@107 | 205 %   to zero values | 
| ivan@107 | 206 | 
| ivan@107 | 207 SMALL.DL(3)=SMALL_init_DL(); | 
| ivan@107 | 208 | 
| ivan@107 | 209 %   Defining fields needed for dictionary learning | 
| ivan@107 | 210 | 
| ivan@107 | 211 SMALL.DL(3).toolbox = 'SPAMS'; | 
| ivan@107 | 212 SMALL.DL(3).name = 'mexTrainDL'; | 
| ivan@107 | 213 | 
| ivan@107 | 214 %   Type 'help mexTrainDL in MATLAB prompt for explanation of parameters. | 
| ivan@107 | 215 | 
| ivan@107 | 216 SMALL.DL(3).param=struct(... | 
| ivan@107 | 217     'D', SMALL.Problem.initdict,... | 
| ivan@107 | 218     'K', SMALL.Problem.p,... | 
| ivan@107 | 219     'lambda', 2,... | 
| ivan@107 | 220     'iter', 200,... | 
| ivan@107 | 221     'mode', 3, ... | 
| ivan@107 | 222     'modeD', 0); | 
| ivan@107 | 223 | 
| ivan@107 | 224 %   Learn the dictionary | 
| ivan@107 | 225 | 
| ivan@107 | 226 SMALL.DL(3) = SMALL_learn(SMALL.Problem, SMALL.DL(3)); | 
| ivan@107 | 227 | 
| ivan@107 | 228 %   Set SMALL.Problem.A dictionary | 
| ivan@107 | 229 %   (backward compatiblity with SPARCO: solver structure communicate | 
| ivan@107 | 230 %   only with Problem structure, ie no direct communication between DL and | 
| ivan@107 | 231 %   solver structures) | 
| ivan@107 | 232 | 
| ivan@107 | 233 SMALL.Problem.A = SMALL.DL(3).D; | 
| ivan@107 | 234 | 
| ivan@107 | 235 %   Setting up reconstruction function | 
| ivan@107 | 236 | 
| ivan@107 | 237 SMALL.Problem.reconstruct = @(x) ImgDenoise_reconstruct(x, SMALL.Problem); | 
| ivan@107 | 238 | 
| ivan@107 | 239 %   Initialising solver structure | 
| ivan@107 | 240 %   Setting solver structure fields (toolbox, name, param, solution, | 
| ivan@107 | 241 %   reconstructed and time) to zero values | 
| ivan@107 | 242 | 
| ivan@107 | 243 SMALL.solver(3)=SMALL_init_solver; | 
| ivan@107 | 244 | 
| ivan@107 | 245 % Defining the parameters needed for image denoising | 
| ivan@107 | 246 | 
| ivan@107 | 247 SMALL.solver(3).toolbox='ompbox'; | 
| ivan@107 | 248 SMALL.solver(3).name='omp2'; | 
| ivan@107 | 249 SMALL.solver(3).param=struct(... | 
| ivan@107 | 250     'epsilon',Edata,... | 
| ivan@107 | 251     'maxatoms', maxatoms); | 
| ivan@107 | 252 | 
| ivan@107 | 253 %   Denoising the image - find the sparse solution in the learned | 
| ivan@107 | 254 %   dictionary for all patches in the image and the end it uses | 
| ivan@107 | 255 %   reconstruction function to reconstruct the patches and put them into a | 
| ivan@107 | 256 %   denoised image | 
| ivan@107 | 257 | 
| ivan@107 | 258 SMALL.solver(3)=SMALL_solve(SMALL.Problem, SMALL.solver(3)); | 
| idamnjanovic@6 | 259 | 
| idamnjanovic@6 | 260 %% | 
| idamnjanovic@6 | 261 % Plot results and save midi files | 
| idamnjanovic@6 | 262 | 
| idamnjanovic@6 | 263 % show results % | 
| idamnjanovic@6 | 264 | 
| idamnjanovic@6 | 265 SMALL_ImgDeNoiseResult(SMALL); |