Clever Geek Handbook

Matlabの周波数フィルタリングに基づいて画像の鮮明度を改善する

はじめに

これまで、数学的な方法の複雑さの速度、コンピューターシステムのリソースの要件などによって特徴付けられる画像の品質を向上させるために、多くのアルゴリズムが開発されてきました。 さらに、最も単純な方法の1つは、周波数領域と空間領域でのフィルタリングに基づく画像処理です。



基本的な概念

周波数領域と空間領域での線形フィルタリングの基礎は、畳み込み定理です。 実際、周波数領域でのフィルター処理の主なアイデアは、特定の方法でF（u、v）を変更するアクセサリーフィルター関数を選択することです。

たとえば、中心関数F（u、v）で乗算され、F（u、v）の高周波成分を減衰させ、低周波成分を比較的変化させずに残し、ローパスフィルターと呼ばれるアクセサリー関数を持つフィルター。

フィルタリングが2次元フーリエ変換に基づいて実装されている場合、画像とその変換は自動的に周期的と見なされます。

周期関数を畳み込む場合、これらの周期が関数の非ゼロ部分の長さに近い位置にある場合、隣接する周期間で干渉効果が発生します。 このような干渉は通常、リターンエラーまたはオーバーラップエラーと呼ばれ、次のように関数にゼロを追加することで抑制できます。

関数f（x、y）とh（x、y）のサイズがABとC Dであると仮定します。fとgにゼロを追加することで、両方とも次元PとQを持つ2つの拡張関数を作成します。

P≥A+ C + 1およびQ≥B+ D + 1

これらの条件を満たすPとQの最小偶数値を計算します。

上記のアルゴリズムを関数として実装します。



  1. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  2. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  3. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  4. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  5. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  6. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  7. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  8. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  9. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  10. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  11. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  12. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  13. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  14. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  15. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  16. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  17. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



  18. function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .



function PQ = paddedsize(AB, CD, PARAM) if nargin==1 PQ=2*AB; elseif nargin == 2&~ischar(CD) PQ = AB+CD ; PQ= 2 *ceil (PQ/2); elseif nargin == 2 m= max (AB); P=2^nextpow2(2*m); PQ=[P, P]; elseif nargin==3 m=max([AB CD]); P=2^nextpow2(2*m); PQ=[P, P]; else error ( ' '); end end * This source code was highlighted with Source Code Highlighter .







フィルターH（u、v）の伝達関数は、F（u、v）の実数部と虚数部で乗算されます。 関数H（u、v）が実数である場合、複素数の実数部と虚数部に同じ数を掛けると位相角は変化しないため、位相方程式からわかるように、積の位相部分は変化しません。 このようなフィルターは、位相シフトがゼロのフィルターと呼ばれます。



  1. 関数 g = dftfilt （ f、H ）
  2. F = fft2 （ f、 サイズ （ H、 1 ） 、 サイズ （ H、 2 ） ） ;
  3. Gi =H。* F;
  4. g = real （ ifft2 （ Gi ） ） ;
  5. g = g （ 1 ： サイズ （ f、 1 ） 、 1 ： サイズ （ f、 2 ） ） ;
  6. 終わり




空間Sobelフィルターに対応する周波数領域でフィルターHを生成します。これにより、画像の垂直エッジが改善されます。



  1. 関数 g = gscale （ f、 varargin ）
  2. 長さ （ varargin ） == 0の場合
  3. method = 'full8' ;
  4. 他に
  5. method = varargin { 1 } ;
  6. 終わり
  7. if strcmp （ class （ f ） 、 'double' ） ＆ （ max （ f （ ： ） ） > 1 | min （ f （ : )） < 0 ）
  8. f = mat2gray （ f ） ;
  9. 終わり
  10. 切り替え方法
  11. ケース 「full8」
  12. g = im2uint8 （ mat2gray （ double （ f ） ） ） ;
  13. ケース 'full16'
  14. g = im2uint16 （ mat2gray （ double （ f ） ） ） ;
  15. ケース 'minmax'
  16. low = varargin { 2 } ; high = varargin { 3 } ;
  17. 低い場合 > 1 | 低< 0 | 高> 1 | 高< 0
  18. エラー （ 「低パラメーターと高パラメーターを変更する必要があります」 ）
  19. 終わり
  20. if strcmp （ クラス （ f ） 、 'double' ）
  21. low_in = min （ f （ :) ） ;
  22. high_in = max （ f （ :) ） ;
  23. elseif strcmp （ クラス （ f ） 、 'uint8' ）
  24. low_in = min （ f （ : )） ./ 255 ;
  25. high_in = max （ f （ : )） ./ 255 ;
  26. elseif strcmp （ クラス （ f ） 、 'uint16' ）
  27. low_in = min （ f （ : )） ./ 65535 ;
  28. high_in = max （ f （ : )） ./ 65535 ;
  29. 終わり
  30. g = imadjust （ f、 [ low_in、high_in ] 、 [ low、high ] ） ;
  31. そうでなければ
  32. エラー （ '無効なメソッド' ）
  33. 終わり




次のdftuv関数は、距離および他の同様のアクションの計算に使用されるグリッド配列を作成します。



  1. 関数 [ U、V ] = dftuv （ M、N ）
  2. u = 0 ： （ M ） ;
  3. v = 0 ： （ N ） ;
  4. idx = find （ u> M / 2 ） ;
  5. u （ idx ） = u （ idx ） ;
  6. idy = find （ v> N / 2 ） ;
  7. v （ idy ） = v （ idy ） ;
  8. [ V、U ] = meshgrid （ v、u ） ;
  9. 終わり




ローパスフィルターを形成する関数を作成する例を示します。



  1. 関数 [ H、D ] = lpfilter （ type 、M、N、D0、n ）
  2. [ U、V ] = dftuv （ M、N ） ;
  3. D = sqrt （ U. ^ 2 + V. ^ 2 ） ;
  4. スイッチ タイプ
  5. ケース 「理想」
  6. H = double （ D <= D0 ） ;
  7. ケース 'btw'
  8. nargin == 4の 場合
  9. nは1です。
  10. 終わり
  11. H = 1 / （ 1 + （ D / D0 ） 。^ （ 2 * n ） ） ;
  12. ケース 「ガウス」
  13. H = exp （ （ D. ^ 2 ） ./ （ 2 * （ D0 ^ 2 ） ） ） ;
  14. そうでなければ
  15. エラー （ 「不明なフィルタータイプ」 ） ;
  16. 終わり
  17. 終わり




次に、ハイパスフィルターを形成する関数の作成例を示します。



  1. 関数 H = hpfilter （ type 、M、N、D0、n ）
  2. [ U、V ] = dftuv （ M、N ） ;
  3. D = sqrt （ U. ^ 2 + V. ^ 2 ） ;
  4. スイッチ タイプ
  5. ケース 「理想」
  6. H = double （ D <= D0 ） ;
  7. ケース 'btw'
  8. nargin == 4の 場合
  9. nは1です。
  10. 終わり
  11. H = 1- （ 1. / （ 1 + （ D / D0 ） 。^ （ 2 * n ） ） ）） ;
  12. ケース 「ガウス」
  13. H = 1- （ exp （ （ D. ^ 2 ） ./ （ 2 * （ D0 ^ 2 ） ） ） ）） ;
  14. そうでなければ
  15. エラー （ 「不明なフィルタータイプ」 ） ;
  16. 終わり
  17. 終わり




また、上記の機能を使用して画像の明瞭さを向上させるスクリプト自体。



  1. img = imread （ 'D：\ Photo \ nature.jpg' ） ;
  2. red = img （ ：、：、 1 ） ;
  3. F = fft2 （ img ） ;
  4. S = fftshift （ log （ 1 + abs （ F ） ） ） ;
  5. S = gscale （ S ） ;
  6. imshow （ img ） 、 図 、imshow （ S ） ;
  7. PQ = paddedsize （ サイズ （赤） ） ;
  8. [ U、V ] = dftuv （ PQ （ 1 ） 、PQ （ 2 ） ） ;
  9. D0 = 0.05 * PQ （ 2 ） ;
  10. F = fft2 （赤、PQ （ 1 ） 、PQ （ 2 ） ） ;
  11. H = exp （ - （ U. ^ 2 + V. ^ 2 ） / （ 2 * （ D0 ^ 2 ） ） ） ;
  12. g = dftfilt （赤、H ） ;
  13. 図 、imshow （ fftshift （ H ） 、 [ ] ） ;
  14. 図 、 メッシュ （ H （ 1 ： 10 ： 500、1 ： 10 ： 500 ） ）
  15. 軸 （ [ 0 50 0 50 0 1 ] ）
  16. カラーマップ （ グレー ）
  17. グリッドオフ
  18. 軸オフ
  19. ビュー （ -25、0 ）
  20. H = fftshift （ hpfilter （ ' gaussian ' 、 500、500、50 ） ） ;
  21. メッシュ （ H （ 1 ： 10 ： 500、1 ： 10 ： 500 ） ）
  22. 軸 （ [ 0 50 0 50 0 1 ] ）
  23. カラーマップ （ [ 0 0 0 ] ）
  24. グリッドオフ
  25. 軸オフ
  26. ビュー （ -163、64 ）
  27. 図 、imshow （ H、 [ ] ） ;
  28. PQ = paddedsize （ サイズ （赤） ） ;
  29. D0 = 0.05。* PQ （ 1 ） ;
  30. H = hpfilter （ 「ガウス」 、PQ （ 1 ） 、PQ （ 2 ） 、D0 ） ;
  31. g = dftfilt （赤、H ） ;
  32. figure 、imshow （ g、 [ 0 255 ] ） ;
  33. PQ = paddedsize （ サイズ （赤） ） ;
  34. D0 = 0.05。* PQ （ 1 ） ;
  35. HBW = hpfilter （ 'btw' 、PQ （ 1 ） 、PQ （ 2 ） 、D0、2 ） ;
  36. H = 0.5 + 2 * HBW;
  37. gbw = dftfilt （赤、HBW ） ;
  38. gbw = gscale （ gbw ） ;
  39. ghf = dftfilt （赤、H ） ;
  40. gbf = gscale （ ghf ） ;
  41. img1 = histeq （ red、 256 ） ;
  42. ghe = histeq （ ghf、 256 ） ;
  43. 図 、imshow （赤） ;
  44. figure 、imshow （ ghe、 [ ] ） ;
  45. figure 、imshow （ img1、 [ ] ） ;




****仕事の例****





































最終結果



More articles:


All Articles