📄 index.html

📁 信号处理系列导航
💻 HTML
字号:
<!DOCTYPE html  PUBLIC "-//W3C//DTD XHTML 1.0 Strict//EN"><html xmlns:mwsh="http://www.mathworks.com/namespace/mcode/v1/syntaxhighlight.dtd">   <head>      <meta http-equiv="Content-Type" content="text/html; charset=utf-8">         <!--This HTML is auto-generated from an M-file.To make changes, update the M-file and republish this document.      -->      <title>Introduction to Wavelet Processing</title>      <meta name="generator" content="MATLAB 7.4">      <meta name="date" content="2008-10-15">      <meta name="m-file" content="index">      <LINK REL="stylesheet" HREF="style.css" TYPE="text/css">   </head>   <body>      <div class="content">         <h1>Introduction to Wavelet Processing</h1>         <introduction>            <p>This numerical tour explore 1D multiresolution analysis.</p>         </introduction>         <h2>Contents</h2>         <div>            <ul>               <li><a href="#1">Installing toolboxes and setting up the path.</a></li>               <li><a href="#8">1D Haar transform</a></li>               <li><a href="#17">Higher order wavelet transform</a></li>            </ul>         </div>         <h2>Installing toolboxes and setting up the path.<a name="1"></a></h2>         <p>You need to download the <a href="../toolbox_general.zip">general purpose toolbox</a> and the <a href="../toolbox_signal.zip">signal toolbox</a>.         </p>         <p>You need to unzip these toolboxes in your working directory, so that you have <tt>toolbox_general/</tt> and <tt>toolbox_signal/</tt> in your directory.         </p>         <p><b>For Scilab user:</b> you must replace the Matlab comment '%' by its Scilab counterpart '//'.         </p>         <p><b>Recommandation:</b> You should create a text file named for instance <tt>numericaltour.sce</tt> (in Scilabe) or <tt>numericaltour.m</tt> to write all the Scilab/Matlab command you want to execute. Then, simply run <tt>exec('numericaltour.sce');</tt> (in Scilab) or <tt>numericaltour;</tt> (in Matlab) to run the commands.         </p>         <p>Execute this line only if you are using Matlab.</p><pre class="codeinput">getd = @(p)path(path,p); <span class="comment">% scilab users must *not* execute this</span></pre><p>Then you can add these toolboxes to the path.</p><pre class="codeinput"><span class="comment">% Add some directories to the path</span>getd(<span class="string">'toolbox_signal/'</span>);getd(<span class="string">'toolbox_general/'</span>);</pre><h2>1D Haar transform<a name="8"></a></h2>         <p>The Haar transform is the simplest orthogonal wavelet transform. It is computed by iterating difference and averaging between            odd and even samples of the signal.         </p>         <p>First we load and display a 1D signal.</p><pre class="codeinput">name = <span class="string">'piece-regular'</span>;n = 512;f = load_signal(name, n);plot(f); axis <span class="string">tight</span>;</pre><img vspace="5" hspace="5" src="index_01.png"> <p>It is easy to split a signal into its odd and even samples, and then compute average and differences.</p><pre class="codeinput">a = f;d = a(1:2:n) - a(2:2:n); <span class="comment">% difference</span>a = a(1:2:n) + a(2:2:n); <span class="comment">% average</span>clf;plot(d); axis(<span class="string">'tight'</span>);</pre><img vspace="5" hspace="5" src="index_02.png"> <p><i>Exercice 1:</i> (the solution is <a href="../private/wavelet_introduction/exo1.m">exo1.m</a>) Implement a full wavelet transform that extract iteratively wavelet coefficients <tt>d</tt> and pack them into a transformed vector <tt>fw</tt> of same size as <tt>f</tt> (you need also to keep the remaining low pass coefficient <tt>a</tt> of size 1). To compute a wavelet transform, you need to iterate this computation on the average vector <tt>a</tt>. Important: for the transform to be orthogonal, you need to multiply the difference and average by <tt>1/sqrt(2)</tt>. Verify that the transform is orthogonal, which means that the energy of the coefficient is the same as the energy of the            signal.         </p><pre class="codeinput">exo1;</pre><pre class="codeoutput">Energy of the signal       = 165562.4476Energy of the coefficients = 165562.4476</pre><img vspace="5" hspace="5" src="index_03.png"> <p>Can you recognize where are the low frequencies and the high frequencies ? You can use the function <tt>plot_wavelet</tt> to help you.         </p><pre class="codeinput">clf; plot_wavelet(fw);</pre><img vspace="5" hspace="5" src="index_04.png"> <p>Inversing the haar transform means retrieving a signal <tt>f1</tt> from the coefficients <tt>fw</tt>. If <tt>fw</tt> are exactely the coefficients of <tt>f</tt>, then <tt>f=f1</tt>.         </p>         <p><i>Exercice 2:</i> (the solution is <a href="../private/wavelet_introduction/exo2.m">exo2.m</a>) Write the inverse wavelet transform that computes <tt>f1</tt> from the coefficients <tt>fw</tt>. Compare <tt>f1</tt> with <tt>f</tt>.         </p><pre class="codeinput">exo2;</pre><pre class="codeoutput">Error |f-f1|/|f| = 8.9882e-16</pre><p>To perform image processing, one usually modify the coefficients <tt>fw</tt>. For instance, it is possible to threshold the coefficients in order to set to 0 the smallest coefficients.         </p><pre class="codeinput">T = 10; <span class="comment">% threshold value</span>fwT = fw .* (abs(fw)&gt;T);clf;subplot(2,1,1);plot_wavelet(fw); axis(<span class="string">'tight'</span>); title(<span class="string">'Original coefficients'</span>);subplot(2,1,2);plot_wavelet(fwT); axis(<span class="string">'tight'</span>); title(<span class="string">'Thresholded coefficients'</span>);</pre><img vspace="5" hspace="5" src="index_05.png"> <p><i>Exercice 3:</i> (the solution is <a href="../private/wavelet_introduction/exo3.m">exo3.m</a>) Find the threshold <tt>T</tt> so that the number of remaining coefficients in <tt>fwT</tt> are <tt>.2*n</tt>. Use this threshold to compute <tt>fWT</tt> and then display the corresponding approximation <tt>f1</tt> of <tt>f</tt>. Try with a smaller number of non zero coefficientS.         </p><pre class="codeinput">exo3;</pre><img vspace="5" hspace="5" src="index_06.png"> <h2>Higher order wavelet transform<a name="17"></a></h2>         <p>Haar approximation give poor results because the Haar wavelet does not have vanishing moments. Furthermore, the reconstructed            image is visually not good looking because the Haar wavelet is not regular. Both issues are solved by using wavelet transform            with a larger number of vanishing moments.         </p>         <p>The function <tt>perform_wavelet_transf</tt> implement a 7/9 biorthogonal wavelet transform. It is close to orthogonal, and has 4 vanishing moments. You can use it to            perform approximation by also thresholding the wavelets coefficients.         </p><pre class="codeinput">Jmin = 2; <span class="comment">% minimum wavelet scale</span>fw = perform_wavelet_transf(f,Jmin,+1);fwT = perform_thresholding(fw, m, <span class="string">'largest'</span>);f2 = perform_wavelet_transf(fwT,Jmin,-1);</pre><p>The comparison with the Haar approximation shows that it gives much better result. This is confirmed by the SNR.</p><pre class="codeinput"><span class="comment">% compute the error</span>ehaar   = snr(f,f1);ebiorth = snr(f,f2);<span class="comment">% display</span>clf;subplot(2,1,1);plot(f1); axis(<span class="string">'tight'</span>); title(strcat([<span class="string">'Haar, SNR='</span> num2str(ehaar) <span class="string">'dB'</span>]));subplot(2,1,2);plot(f2); axis(<span class="string">'tight'</span>); title(strcat([<span class="string">'7/9, SNR='</span> num2str(ebiorth) <span class="string">'dB'</span>]));</pre><img vspace="5" hspace="5" src="index_07.png"> <p class="footer"><br>            Copyright  &reg; 2008 Gabriel Peyre<br></p>      </div>      <!--##### SOURCE BEGIN #####%% Introduction to Wavelet Processing% This numerical tour explore 1D multiresolution analysis.%% Installing toolboxes and setting up the path.%%% You need to download the % <../toolbox_general.zip general purpose toolbox>% and the <../toolbox_signal.zip signal toolbox>.%%% You need to unzip these toolboxes in your working directory, so% that you have |toolbox_general/| and |toolbox_signal/| in your% directory.%%% *For Scilab user:* you must replace the Matlab comment '%' by its Scilab% counterpart '//'.%%% *Recommandation:* You should create a text file named for instance% |numericaltour.sce| (in Scilabe) or |numericaltour.m| to write all the% Scilab/Matlab command you want to execute. Then, simply run% |exec('numericaltour.sce');| (in Scilab) or |numericaltour;| (in Matlab)% to run the commands. %%% Execute this line only if you are using Matlab.getd = @(p)path(path,p); % scilab users must *not* execute this%%% Then you can add these toolboxes to the path.% Add some directories to the pathgetd('toolbox_signal/');getd('toolbox_general/');%% 1D Haar transform% The Haar transform is the simplest orthogonal wavelet transform. It is% computed by iterating difference and averaging between odd and even% samples of the signal.%%% First we load and display a 1D signal.name = 'piece-regular';n = 512;f = load_signal(name, n);plot(f); axis tight;%%% It is easy to split a signal into its odd and even samples, and then% compute average and differences.a = f;d = a(1:2:n) - a(2:2:n); % difference a = a(1:2:n) + a(2:2:n); % averageclf;plot(d); axis('tight');%%% _Exercice 1:_ (the solution is <../private/wavelet_introduction/exo1.m exo1.m>)% Implement a full wavelet transform that extract iteratively wavelet% coefficients |d| and pack them into a transformed vector |fw| of same% size as |f| (you need also to keep the remaining low pass coefficient% |a| of size 1).% To compute a wavelet transform, you need to iterate this computation on% the average vector |a|.% Important: for the transform to be orthogonal, you need to multiply the% difference and average by |1/sqrt(2)|. Verify that the transform is% orthogonal, which means that the energy of the coefficient is the same% as the energy of the signal.exo1;%%% Can you recognize where are the low frequencies and the high frequencies% ? You can use the function |plot_wavelet| to help you.clf; plot_wavelet(fw);%%% Inversing the haar transform means retrieving a signal |f1| from the% coefficients |fw|. If |fw| are exactely the coefficients of |f|, then% |f=f1|.  %%% _Exercice 2:_ (the solution is <../private/wavelet_introduction/exo2.m exo2.m>)% Write the inverse wavelet transform that computes |f1| from the% coefficients |fw|. Compare |f1| with |f|. exo2;%%% To perform image processing, one usually modify the coefficients |fw|.% For instance, it is possible to threshold the coefficients in order to% set to 0 the smallest coefficients. T = 10; % threshold valuefwT = fw .* (abs(fw)>T);clf;subplot(2,1,1);plot_wavelet(fw); axis('tight'); title('Original coefficients');subplot(2,1,2);plot_wavelet(fwT); axis('tight'); title('Thresholded coefficients');%%% _Exercice 3:_ (the solution is <../private/wavelet_introduction/exo3.m exo3.m>)% Find the threshold |T| so that the number of remaining coefficients in% |fwT| are |.2*n|. Use this threshold to compute |fWT| and then display% the corresponding approximation |f1| of |f|. Try with a smaller number% of non zero coefficientS.exo3;%% Higher order wavelet transform% Haar approximation give poor results because the Haar wavelet does not% have vanishing moments. Furthermore, the reconstructed image is visually% not good looking because the Haar wavelet is not regular. Both issues are% solved by using wavelet transform with a larger number of vanishing% moments.%%% The function |perform_wavelet_transf| implement a 7/9 biorthogonal% wavelet transform. It is close to orthogonal, and has 4 vanishing% moments. You can use it to perform approximation by also thresholding the% wavelets coefficients.Jmin = 2; % minimum wavelet scalefw = perform_wavelet_transf(f,Jmin,+1);fwT = perform_thresholding(fw, m, 'largest');f2 = perform_wavelet_transf(fwT,Jmin,-1);%%% The comparison with the Haar approximation shows that it gives much% better result. This is confirmed by the SNR.% compute the errorehaar   = snr(f,f1);ebiorth = snr(f,f2);% displayclf;subplot(2,1,1);plot(f1); axis('tight'); title(strcat(['Haar, SNR=' num2str(ehaar) 'dB']));subplot(2,1,2);plot(f2); axis('tight'); title(strcat(['7/9, SNR=' num2str(ebiorth) 'dB']));##### SOURCE END #####-->   </body></html>
⌨️ 快捷键说明

复制代码 Ctrl + C
搜索代码 Ctrl + F
全屏模式 F11
切换主题 Ctrl + Shift + D
显示快捷键 ?
增大字号 Ctrl + =
减小字号 Ctrl + -