画像情報(輝度データ)の変換


平滑化による画像データの変換

画像データの平滑化とはあるピクセルに注目しその輝度をそのピクセルとその周りをぐるりと囲むピクセルの輝度の平均、あるいは重み付き平均をとって平均化させようというものである。これによりさまざまな雑音を低減させようとするものである。いわゆる画像フィルタの1種である。今回この教材で使用する画像データは以下のようなものである。これをさっそく自分の作業ディレクトリにコピーしておこう。

> cp ~fukunaga/ipwork/tmuinthouse2.pgm ./

ノイズの付加(シミュレーション)

ディジタル画像上の加工過程で現れる雑音(ノイズ)はインパルスノイズとガウシアンノイズと呼ばれるものに分類される事が知られています。この2つのノイズを人工的に作り出し画像を劣化させる試みをします。

インパルスノイズ

位置に独立でランダムにある頻度をもって現れるノイズです。あるピクセルでノイズが載るとそのピクセルデータはゼロ(黒とび)になるか最大値(白とび)になってしまいます。

適当な頻度でインパルスノイズを付加するプログラムを以下に示します。確率pを定義し、その確率のさらに1/2でそのピクセルの現在値には関係なくそのピクセルの輝度を0、あるいはMax Pixcel Valueととるようにします。

実行時の初期パラメータを2つとり、最初がオリジナル画像ppmファイル名、2番目がノイズレートです。プログラムの下に載せた図はオリジナルに対して5%のノイズを(さらにそのうち50%が黒点、50%が白点)発生させて劣化させたものです。このプログラムをよく読み実行させてください。

> cp ~fukunaga/ipwork/ppmImpuseAdd.java ./

> javac ppmImpulseAdd.java

//import ppmread.class ;
class ppmImpuseAdd {
    public static void main(String args[]) {
        int i,j,k ;
        int histdata[] ;
        double nratio, fwork ; 
        ppmread image = new ppmread() ;
        if(args.length < 2) {
            System.out.println("Usage: java ppmImpuseAdd filename Noise_ratio");
            System.exit(1) ;
        }
        image.ppmmain(args[0]) ;
        nratio = Double.parseDouble(args[1]) ;
        if(nratio<0.0 || nratio>=1.0) {
            System.out.println("Noise Ratio should be in between 0.0 and 1.0");
            System.exit(1) ;
        }
        // PGM Header output 
        System.out.println("P"+image.magicidno) ;
        System.out.println("# "+image.sfilename+" is modified Brighter") ;
        System.out.print(image.width+" "+image.height) ;
        if(image.magicidno >1 ) System.out.println(" "+image.maxval) ;
        System.out.println("# Data created by ppmHistdata") ;
        System.out.println("# input filename = "+ image.sfilename ) ;
        System.out.println("# Width    = "+image.width) ;
        System.out.println("# Height   = "+image.height) ;
        if(image.magicidno>1) System.out.println("# Maxval   = "+image.maxval) ;
        
        for(i=0;i<image.height;i++) 
            for(j=0;j<image.width;j++) {
                if(Math.random()<nratio) 
                  if(Math.random()>0.5) image.gdata[j][i] = 0 ;
                  else image.gdata[j][i] = (short)image.maxval ;
            }
        int ndataMaxLine=50 ;
        int rdataline,rdatalast ;
        rdataline = image.width/ndataMaxLine ;
        rdatalast = image.width%ndataMaxLine ;
        for(k=0;k<image.height;k++) {
            for(i=0;i<rdataline;i++) {
                for(j=0;j<ndataMaxLine;j++) {
                    int ij=i*ndataMaxLine+j ;
                    System.out.print(image.gdata[ij][k]+" ") ;

                }
                System.out.println() ;
            }
            for(j=0;j<rdatalast;j++) {
                int ij=rdataline*ndataMaxLine+j ;
                System.out.print(image.gdata[ij][k]+" ") ;
            }
            System.out.println() ;
        }
    }
}

> java ppmImpulseAdd tmuinthouse2.pgm 0.05 > noisetmuinthouse2.pgm

この操作でできた画像が以下のもの。

ガウシアンノイズ

各ピクセルごとに現在の輝度を中心に適当な標準偏差σ(20〜40)で正規分布乱数をふり輝度を変更させます。正規分布の特性により輝度変化が現在値±σ内に収まる確率が68%、あるいは現在値±2σより大となる確率は4.55%になります。したがって480×360ピクセルの画像では7860余りののセルが大きな輝度変化を起こしてしまう可能性があります。以下のプログラム第2パラメータに実数で標準偏差を指定します。

> cp ~fukunaga/ipwork/ppmGaussianAdd.java ./

> javac ppmGaussianAdd.java

//import ppmread.class ;
class ppmGaussianAdd {
    public static double gaussRandom() {
        //  Gaussian Random Number generator (mean =0 ,sigma=1)
        int i ;
        double x ;
        x=0.0 ;
        for(i=0;i<12;i++) x += Math.random() ;
        x -=6.0 ;
        return(x) ;
    }
        
    public static void main(String args[]) {
        int i,j,k ;
        int histdata[] ;
        double nratio, fwork ; 
        double sigma ;
        double rtwopi1=0.39894228 ;
        double pgaussmax ;
        ppmread image = new ppmread() ;
        if(args.length < 2) {
            System.out.println("Usage: java ppmGaussianAdd filename SIGMA");
            System.exit(1) ;
        }
        image.ppmmain(args[0]) ;
        sigma = Double.parseDouble(args[1]) ;
        if(sigma<0.0 ) {
            System.out.println("sigma must be > 0");
            System.exit(1) ;
        }
        // PGM Header output 
        System.out.println("P"+image.magicidno) ;
        System.out.println("# "+image.sfilename+" is modified Brighter") ;
        System.out.print(image.width+" "+image.height) ;
        if(image.magicidno >1 ) System.out.println(" "+image.maxval) ;
        System.out.println("# Data created by ppmGaussianAdd") ;
        System.out.println("# input filename = "+ image.sfilename ) ;
        System.out.println("# Width    = "+image.width) ;
        System.out.println("# Height   = "+image.height) ;
        if(image.magicidno>1) System.out.println("# Maxval   = "+image.maxval) ;
        
        for(i=0;i<image.height;i++) 
            for(j=0;j<image.width;j++) {
                image.gdata[j][i] += gaussRandom()*sigma ;
                if(image.gdata[j][i]<0) image.gdata[j][i]=0 ;
                if(image.gdata[j][i]>image.maxval) 
                    image.gdata[j][i]=(short)image.maxval ;
            }
        int ndataMaxLine=50 ;
        int rdataline,rdatalast ;
        rdataline = image.width/ndataMaxLine ;
        rdatalast = image.width%ndataMaxLine ;
        for(k=0;k<image.height;k++) {
            for(i=0;i<rdataline;i++) {
                for(j=0;j<ndataMaxLine;j++) {
                    int ij=i*ndataMaxLine+j ;
                    System.out.print(image.gdata[ij][k]+" ") ;
                }
                System.out.println() ;
            }
            for(j=0;j<rdatalast;j++) {
                int ij=rdataline*ndataMaxLine+j ;
                System.out.print(image.gdata[ij][k]+" ") ;
            }
            System.out.println() ;
        }
    }
}

> java ppmGaussianAdd tmuinthouse2.pgm 20.0 > gnoisetmuinthouse2.pgm

この操作で得られた画像は以下のようです。

単純平均化

あるピクセル(i,j)の元データをf(i,j)、変換後のデータをg(i,j)とすると(i,j)ピクセルおよびその周りのピクセル値をとり3×3行列を作る。そして以下の図のごとくその行列に重み行列をかけてやることにより自由自在のフィルタを作成させてやることができる。

そしてこの操作をエッジをのぞきすべてに適用する。一般的にこのような画像の局所ピクセルとその大きさに対応する行列をかけて画像データの変換をすることを局所積和(平均)と呼びます。画像データ変換で常に使われるテクニックです。

局所平均フィルタ

以下の行列を使った場合、もっとも単純な平滑化の手法として広く使われている。しかしノイズ(雑音)のぼかしとともに濃淡画像をも滑らかにする。輪郭もあいまいになる。

局所加重平均フィルタ

重み係数行列としてよく使われるのが以下のようなものである。どちらも中心画素に近いほどその影響を高めるように配慮している。

以下のプログラムは上記の3つのマトリックスのどれかを実行時に指定することにより(第2番目の初期パラメータを0,1,あるいは2)として3つの平滑化プロセスを1つのプログラムで実現できるようにしたものである。第1番目のパラメータは対象となる画像pgmファイルである。

> cp ~fukunaga/ipwork/ppmImagSmooth.java ./

> javac ppmImagSmooth.java

//import ppmread.class ;
class ppmImagSmooth {
    public static void main(String args[]) {
        int i,j,k ;
        short imageData[][] ;
        ppmread objImage = new ppmread() ;
        double mfactor ;
        short weightMethod=0 ;
        double W[][] = new double[3][3] ;
        double f[][] = new double[3][3] ;
        double Wval ;
        double[][] W0 ={{1.0,1.0,1.0},
                        {1.0,1.0,1.0},
                        {1.0,1.0,1.0}} ;
        double W0val = 9.0 ;
        double[][] W1 ={{0.1,0.1,0.1},
                        {0.1,0.2,0.1},
                        {0.1,0.1,0.1}} ;
        double W1val = 1.0 ;
        double[][] W2 ={{0.0625,0.125,0.0625},
                        {0.125,0.25,0.125},
                        {0.0625,0.125,0.0625}} ;
        double W2val = 1.0 ;
        if(args.length<1) {
            System.out.println
                ("Usage: java ppmImagSmooth filename weight_method") ;
            System.exit(0) ;
        }
        objImage.ppmmain(args[0]) ;
        if(args.length > 1)
             weightMethod  = (short)Integer.parseInt(args[1]) ;
        if(weightMethod < 0 || weightMethod > 2) {
            System.out.println
                ("weight_method must be in betwee 0 and 2") ;
            System.exit(0) ;
        }
        Wval = 9.0 ;
        switch(weightMethod) {
        case 0:
            Wval = W0val ;
            break ;
        case 1:
            Wval = W1val ;
            break ;
        case 2:

            Wval = W2val ;
            break;
        }
        for(i=0;i<3;i++) for(j=0;j<3;j++) 
            switch(weightMethod) {
            case 0:
                W[i][j]=W0[i][j] ;
                break ;
            case 1:
                W[i][j]=W1[i][j] ;
                break ;
            case 2:
                W[i][j]=W2[i][j] ;
                break;
            }
        imageData=new short [objImage.width][objImage.height] ;
        System.out.println("P"+objImage.magicidno) ;
        System.out.println("# Data is created by ppmImagSmooth") ;
        System.out.println("# Weight Method selected is "+weightMethod) ;
        System.out.println("# filename = "+ objImage.sfilename ) ;
        System.out.println("# Width    = "+objImage.width) ;
        System.out.println("# Height   = "+objImage.height) ;
        if(objImage.magicidno>1) 
            System.out.println("# Maxval   = "+objImage.maxval) ;
        System.out.print(objImage.width+" "+objImage.height) ;
        if(objImage.magicidno >1 ) System.out.println(" "+objImage.maxval) ;
        double work ;
        int ii,jj ;
        for(i=1;i<objImage.height-1;i++) 
            for(j=1;j<objImage.width-1;j++) {
                for(ii=-1;ii<2;ii++) 
                    for(jj=-1;jj<2;jj++)
                        f[ii+1][jj+1]=objImage.gdata[j+jj][i+ii] ;
                work =0.0 ;
                for(ii=0;ii<3;ii++) 
                    for(jj=0;jj<3;jj++) 
                        work += W[ii][jj]*f[jj][ii] ;
                work /= Wval ;
                imageData[j][i] = (short)Math.rint(work) ;
            }
        int ndataMaxLine=50 ;
        int rdataline,rdatalast ;
        rdataline = objImage.width/ndataMaxLine ;
        rdatalast = objImage.width%ndataMaxLine ;

        for(k=0;k<objImage.height;k++) {
            for(i=0;i<rdataline;i++) {
                for(j=0;j<ndataMaxLine;j++) {
                    int ij=i*ndataMaxLine+j ;
                    System.out.print(imageData[ij][k]+" ") ;
                }
                System.out.println() ;
            }
            for(j=0;j<rdatalast;j++) {
                int ij=rdataline*ndataMaxLine+j ;
                System.out.print(imageData[ij][k]+" ") ;
            }
            System.out.println() ;
        }
    }
}

このプログラムの実行結果を以下に示す。

>java ppmImagSmooth tmuinthouse2.pgm 0

>java ppmImagSmooth noisetmuinthouse2.pgm 0

>java ppmImagSmooth noisetmuinthouse2.pgm 2

>java ppmImagSmooth gnoisetmuinthouse2.pgm 0

画質をそこなわない平滑化

メディアンフィルタ

局所(加重)平均フィルタはノイズとともにエッジもぼかしてしまうが、このメディアンフィルタではエッジのぼかしは押さえられ、ノイズのみをぼかしてしまうすぐれたものです。局所領域のメディアン値をその局所ピクセルの値としてしまうものです。局所領域(3×3)で一瞬にしてヒストグラムを作らなければなりません。

メディアンはある分布しているデータをエントリ数が等しくなるように等分し、その等分割する縦軸のx(横)軸と交わる点のx値をとったものです。つまりヒストグラムを作ってそのグラフ中の全エントリ数の1/2になる横軸をとるものです。また言い方を変えればそれはデータを昇順に配列(ソート)し、その中位を占める値のことです。上の例で行くと9個のエントリの並び、

5 3 4 3 10 5 3 4 5

を昇順ソートし

3 3 3 4 4 5 5 5 10

とした場合の5番目のデータ4をとるということです。以下のプログラムではこのソートにクイックソートアルゴリズムを利用して膨大なデータ量(局所範囲を3×3に限定しても、480×360の画像だと9個のソートを17万回以上を行わなければならない)処理を高速化しています。

以下にそのメディアンフィルタのプログラムとその適用結果を示します。

なおこのプログラムではclass quickSortを使います。quickSort.classファイルが同じ作業ディレクトリになければなりません。以下の要領でこのクラスを作成しておいてください(ここにはそのソースを掲げませんがプログラムの内容は理解しておいてください。再帰処理を使っています)。

> cp ~fukunaga/ipwork/quickSort.java ./

> javac quickSort.java

これらの作業が終了した後(なんの警告、注意メッセージも出なかったことを確認して)、さらに

> cp ~fukunaga/ipwork/ppmImagMedianSmooth.java ./

> javac ppmMedianSmooth.java

//import ppmread.class ;
//import quickSort.class ;
class ppmImagMedianSmooth {
    public static void main(String args[]) {
        int i,j,k,n ;
        int ii,jj ;
        short imageData[][] ;
        int []tobesorted = new int[9] ;
        ppmread objImage = new ppmread() ;
        if(args.length<1) {
            System.out.println
                ("Usage: java ppmImagMedianSmooth filename") ;
            System.exit(0) ;
        }
        quickSort q = new quickSort() ; 
        objImage.ppmmain(args[0]) ;
        imageData=new short [objImage.width][objImage.height] ;
        System.out.println("P"+objImage.magicidno) ;
        System.out.println("# Data is created by ppmImagMedianSmooth") ;
        System.out.println("# filename = "+ objImage.sfilename ) ;
        System.out.println("# Width    = "+objImage.width) ;
        System.out.println("# Height   = "+objImage.height) ;
        if(objImage.magicidno>1) 
            System.out.println("# Maxval   = "+objImage.maxval) ;
        System.out.print(objImage.width+" "+objImage.height) ;
        if(objImage.magicidno >1 ) System.out.println(" "+objImage.maxval) ;
        for(i=1;i<objImage.height-1;i++) 
            for(j=1;j<objImage.width-1;j++) {
                n = 0 ;
                for(ii=-1;ii<2;ii++) 
                    for(jj=-1;jj<2;jj++)
                        tobesorted[n++]=objImage.gdata[j+jj][i+ii] ;
                //sort local image data
                q.quicksort(tobesorted,0,8) ;            
                //center of the local is replaced with median value 
                imageData[j][i] = (short)tobesorted[5] ;  
            }
        int ndataMaxLine=50 ;
        int rdataline,rdatalast ;
        rdataline = objImage.width/ndataMaxLine ;
        rdatalast = objImage.width%ndataMaxLine ;
        for(k=0;k<objImage.height;k++) {
            for(i=0;i<rdataline;i++) {
                for(j=0;j<ndataMaxLine;j++) {

                    int ij=i*ndataMaxLine+j ;
                    System.out.print(imageData[ij][k]+" ") ;
                }
                System.out.println() ;
            }
            for(j=0;j<rdatalast;j++) {
                int ij=rdataline*ndataMaxLine+j ;
                System.out.print(imageData[ij][k]+" ") ;
            }
            System.out.println() ;
        }
    }
}

そして例えば以下のように実行します。

> java ppmMedianSmooth tmuinthouse2.pgm > medtmuinthouse2.pgm

いくつか実行結果を示します。

>java ppmImagMedianSmooth tmuinthouse2.pgm (オリジナル画像への適用)

>java ppmImagMedianSmooth imptmuinthouse2.pgm (インパルスノイズを付加した画像への適用)

実習:

ガウシアンノイズを付加した画像に対してメディアンフイルタを適用したらどうなるか確認してください。

研究(課題1):

局所平均フィルタのところの下線部を実証するデータを作り出すプログラムを作成してください.それをそのデータから示してください.

研究(課題2):

各種平滑化フィルタによって画像が改善されたことが目でも確認されているが、これを定量的に評価するため平均二乗誤差(Mean Square Error; MSE)による数値検査を導入する:原画の輝度値をf(i,j)、フィルタ適用後のそれをg(i,j)とする.m×nの画素を持つ画像でMSEは

とする.また正規化MSEを以下のように定義する;

これらを求めるプログラムを以下を参考に作ってください.ppgImagMSE.java

//import ppmread.class ;
class ppmImagMSE {
    public static void main(String args[]) {
        int i,j,k ;
        int difffg,sumf ; 
        double MSE,nMSE ;
        ppmread origImage = new ppmread() ;
        ppmread modifImage = new ppmread() ;
        if(args.length < 1) {
            System.out.println
                ("Usage: java ppmImagMSE origFilename modifFilename");
            System.exit(1) ;
        }
        origImage.ppmmain(args[0]) ;
        modifImage.ppmmain(args[1]) ;
        if(origImage.width == modifImage.width)
            if(origImage.height != modifImage.height) {
                System.out.println("different image heights!") ;
                System.exit(1) ;
            }
            else ;
        else {
          System.out.println("different image widths!") ;
          System.exit(1) ;
      }
      //   以下 計算部  //

実行は例えば

> java ppmImagMSE tmuinthouse2.pgm Imphouse2.pgm

などとして

MSE = 178.5518 …
nMSE= 1.1781 …

など標準出力に結果を示す.


研究(課題3):

見た目ではインパルスノイズにはメディアンフィルタ、ガウシアンノイズにはいわゆる局所(加重)平均フィルタがすぐれているようである。とにかくインパルスノイズに平均フィルタは使えない。なぜインパルス→メディアン、ガウシアン→平均フィルタなのかそのノイズの特性、平滑化アルゴリズムの特性を考慮して理由を述べてください.

研究(課題4):

関西大学の福西潤一氏は方向性差分フィルタを提唱しインパルスノイズ除去に効果があると報告している
http://joho.densi.kansai-u.ac.jp/isp/impulse-j.html).以下のようなアルゴリズムである.まず3×3のピクセルで方向を以下の4つと定めておく.

1)各方向の差分の計算

2)各方向の平均値

3)ピクセル(i,j)の濃度値を上式のgijとする.つまり差分がもっとも小さくなるnを選びそのm[n]を更新値とするのである.

このアルゴリズムを実装せよ.プログラムの名前はppmImagSmoothDD.javaとする.これをインパルスノイズを加えた画像(tmuinthouse2.pgm+impulse noise 10%)に(1回)適用し結果をppmImagSmoothと比較せよ.定量的な評価を研究(課題2)のプログラムppmImagMSEを用いて行え.


エッジ・線の検出

エッジ要素の抽出

平滑化のところで見た局所(3×3行列)積和計算技法を用いて画像の線、輪郭を引き出すことができる。一般にエッジと呼ばれるところは明るさが急激的に変化しているところで、その明るさの変化をイラストしてみると以下の図のように画像上でだいたい3種類の構造となって確認されます。

したがってこのグラディエントあるいはラプラシアンを計算してやればエッジの検出が原理的に可能です。グラディエントは1次微分、ラプラシアンは2次微分をとるという事に対応します。

x方向、y方向の単位ベクトルをux、uyとするとグラディエント(2次元)はベクトルであたえられ

この強度(微分値)を以下のように定義し、実際の変換後のピクセル(i,j)のデータ値と一般的にしているようです(変形も考案されている)。

ラプラシアンは

となります。

画像データのようなグリッド(格子)点での微分は差分で計算されます。x、y方向の1次、2次差分は(現在ピクセル(i,j)点にいるとして)以下のようにあたえられます。1次微分をとるとき自分自身の値は使わず隣接する前後、上下のピクセルのデータを用いることにします。

fx = (-1) f(i-1,j) + (0) f(i,j) + (+1)f(i+1,j),

より係数項行ベクトル [-1,0,1]

fy = (-1) f(i,j-1) + (0) f(i,j) + (+1)f(i,j+1)

より係数項列ベクトル(転置) [ -1,0,1]

を書き出して、さらに3×3ピクセルに拡張させるとそれぞれの重み行列は以下のようになります。

上のfx, fyも今までの重み行列の記法とあわせて、Wと呼んでいくことにする。

同じように2次差分係数も以下の式より

係数項行(列)ベクトル=[1 -2 1]となるので fxx+fyy でラプラシアンとすれば、その重み行列は以下のようになります。つまりfxxおよびfyyは以下のように与えられるので、

あるいは斜め方向(クロス項)も考慮したが利用される。

さて3×3の局所積和平均に関してこの操作は先に示したプログラムppmImagSmooth.javaとまったく同じです。行列Wの数値が違うだけです。このプログラムを拡張すればいいのですが、ここでは一応別プログラムを用意しました。ppmEgdeDetect.javaです。

~fukunaga/ipwork/ppmEdgeDetect.java

//import ppmread.class ;
class ppmEdgeDetect {
    public static void main(String args[]) {
        int i,j,k ;
        short imageData[][] ;
        ppmread objImage = new ppmread() ;
        double mfactor ;
        short weightMethod=0 ;
        double Wx[][] = new double[3][3] ;
        double Wy[][] = new double[3][3] ;
        double f[][] = new double[3][3] ;
        double[][] Wgradx ={{-1,0,1},
                            {-1,0,1},
                            {-1,0,1}} ;
        double[][] Wgrady ={{-1,-1,-1},
                            {0,0,0},
                        {1,1,1}} ;
        double[][] Wsobelx ={{-1,0,1},
                             {-2,0,2},
                             {-1,0,1}} ;
        double[][] Wsobely ={{-1,-2,-1},
                             {0,0,0},
                             {1,2,1}} ;
        double[][] Wlap ={{0,1,0},
                          {1,-4,1},
                        {0,1,0}} ;
        double[][] Wlapcros ={{1,1,1},
                             {1,-8,1},
                             {1,1,1}} ;
        if(args.length<1) {
            System.err.println
                ("Usage: java ppmEdgeDetect filename weight_method") ;
            System.exit(0) ;
        }
        if(args.length > 1)
             weightMethod  = (short)Integer.parseInt(args[1]) ;
        if(weightMethod < 0 || weightMethod > 7) {
            System.err.println
                ("\tNo.\tWeighted Matrix") ;
            System.err.println
                ("\t0\tGradient x y" ) ;
            System.err.println
                ("\t1\tGradient x") ;
            System.err.println
                ("\t2\tGradient y") ;
            System.err.println
                ("\t3\tSobel x") ;
            System.err.println
                ("\t4\tSobel y") ;
            System.err.println
                ("\t5\tSobel x y") ;
            System.err.println
                ("\t6\tLaplacian") ;
            System.err.println
                ("\t7\tLaplacian (incl. cross term)") ;
            System.exit(1) ;
        }
        if(args.length<1) {
            System.err.println
                ("Usage: java ppmEdgeDetect filename weight_method") ;
            System.exit(0) ;
        }

        objImage.ppmmain(args[0]) ;
        for(i=0;i<3;i++) for(j=0;j<3;j++) 
            switch(weightMethod) {
            case 0: 
                Wx[i][j]=Wgradx[i][j] ;
                Wy[i][j]=Wgrady[i][j] ;
                break ;
            case 1:
                Wx[i][j]=Wgradx[i][j] ;
                break ;
            case 2:
                Wx[i][j]=Wgrady[i][j] ;
                break;
            case 3:
                Wx[i][j]=Wsobelx[i][j] ;
                break;
            case 4:
                Wx[i][j]=Wsobely[i][j] ;
                break;
            case 5:
                Wx[i][j]=Wsobelx[i][j] ;
                Wy[i][j]=Wsobely[i][j] ;
                break;
            case 6:
                Wx[i][j]=Wlap[i][j] ;
                break;
            case 7:
                Wx[i][j]=Wlapcros[i][j] ;
                break;
            }
        //if(weightMethod==0) {
        //  System.err.println("Wy00,Wy11="+Wy[0][0]+" "+Wy[1][1]) ;
        //}
        imageData=new short [objImage.width][objImage.height] ;
        System.out.println("P"+objImage.magicidno) ;
        System.out.println("# Data is created by ppmEdgeDetect") ;
        System.out.println("# Weight Method selected is "+weightMethod) ;
        System.out.println("# filename = "+ objImage.sfilename ) ;
        System.out.println("# Width    = "+objImage.width) ;
        System.out.println("# Height   = "+objImage.height) ;
        if(objImage.magicidno>1) 
            System.out.println("# Maxval   = "+objImage.maxval) ;
        System.out.print(objImage.width+" "+objImage.height) ;
        if(objImage.magicidno >1 ) System.out.println(" "+objImage.maxval) ;
        double workx,worky ;
        int ii,jj,maxdata=0,mindata=0 ;
        for(i=1;i<objImage.height-1;i++) {
            for(j=1;j<objImage.width-1;j++) {
                for(ii=-1;ii<2;ii++) 
                    for(jj=-1;jj<2;jj++)
                        f[ii+1][jj+1]=objImage.gdata[j+jj][i+ii] ;
                workx =0.0 ;
                worky =0.0 ;
                for(ii=0;ii<3;ii++) 
                    for(jj=0;jj<3;jj++) 
                        workx += Wx[ii][jj]*f[jj][ii] ;
                if(weightMethod==0||weightMethod==5) 
                    for(ii=0;ii<3;ii++) 
                        for(jj=0;jj<3;jj++) 
                            worky += Wy[ii][jj]*f[jj][ii] ;
                
                imageData[j][i] = (short)Math.rint(workx) ;
                if(weightMethod==0||weightMethod==5) 
                    imageData[j][i]=
                        (short)Math.rint(Math.sqrt(workx*workx+worky*worky)) ;
                if(imageData[j][i]>255) imageData[j][i]=255;
                if(imageData[j][i]<0) imageData[j][i]=0 ;
                if(maxdata<imageData[j][i]) maxdata=imageData[j][i] ;
                if(mindata>imageData[j][i]) mindata=imageData[j][i] ;
            }
        }
        //      System.err.println("Max = "+maxdata+"\nMin = "+mindata) ;
        int ndataMaxLine=50 ;
        int rdataline,rdatalast ;
        rdataline = objImage.width/ndataMaxLine ;
        rdatalast = objImage.width%ndataMaxLine ;
        for(k=0;k<objImage.height;k++) {
            for(i=0;i<rdataline;i++) {
                for(j=0;j<ndataMaxLine;j++) {
                    int ij=i*ndataMaxLine+j ;
                    System.out.print(imageData[ij][k]+" ") ;
                }
                System.out.println() ;
            }
            for(j=0;j<rdatalast;j++) {
                int ij=rdataline*ndataMaxLine+j ;
                System.out.print(imageData[ij][k]+" ") ;
            }
            System.out.println() ;
        }
    }
}

> cp ~fukunaga/ipwork/ppmEdgeDetect.java ./

> javac ppmEdgeDetect.java

実行は初期パラメータとして入力画像データファイル名、行列の種類を整数であたえます。行列の種類の整数を省略すると0が仮定されます。数字と種類は以下のようです。この整数になんらかの負の値をあたえると正の整数値0から7までの意味合いを表示してくれます。つまり以下の表が出力されます。

入力指定数値 内容 行列
0 グラディエント sqrt(fx2+fy2)
1 グラディエント fx
2 グラディエント fy
3 グラディエント fx (Sobel)
4 グラディエント fy (Sobel)
5 Sobel sqrt(fx2+fy2)
6 ラプラシアン
7 ラプラシアン 含クロス項

tmuinthouse2.pgmにSobelグラディエント行列でxy両方向のエッジ検出を行おうとするなら、

> java ppmEdgeDetect tmuinthouse2.pgm 5 > sobelxytmuinthouse2.pgm

などとすればよいです。エッジ検出の効果を確認してください。

いくつかサンプルを掲げます。

グラディエント x

グラディエント y

グラディエント xy

ラプラシアン(含クロス項)

画像データの幾何学変換

拡大・縮小

拡大・縮小は以下の式で処理される。

つまりu=ax,v=byでaはx方向の拡大(縮小)率、bはy方向のそれである.

以下のプログラムは拡大・縮小を行うプログラムppmScaling.javaである.

課題2

プログラム実行時にa,bをあたえ(ともに実数で、パラメータ1つの場合自動的にx、yともに同じ比率のスケーリングで)原画像を拡大縮小するプログラムppmScaling.javaを作ってください。この際ファイルはグレイスケール(pgm)、カラー(ppm)である場合も考慮してください。

例えばx方向2倍の拡大によってピクセルが増えた場合、増えたピクセルの輝度あるいが色度は決定されていません。

最初

1 2 1
1 2 1
1 2 1

拡大後(新規画素の輝度情報なし)

1 2 1
1 2 1
1 2 1

縮小後(なくなった画素情報の分配なし)

1 2
1 2
1 2

そこで決めるべき輝度データの上下、左右、斜めの位置からの値の平均値(平均値の重みは距離の逆比例)をとって画素の値を再分配するようにします。

そのために以下の線形補間法を利用しました.

上の図より新たに(拡大・縮小)設定された点g(u,v)は従来の画素濃度からの線形補間:

で与えられる.

プログラムを見てみよう.補間はx方向i=0、i=width-1、j=0、j=height-1で例外的取り扱いをする(特に拡大の場合).そこを注意してちょっと長いがプログラムを読み取ってほしい.

//import ppmread.class ;
class ppmScaling {
    public static void main(String args[]) {
        short imageData[][] ;
        int enlargex,enlargey ;
        //int newcoord[][] ;
        int newwidth, newheight ;
        //      byte bar[]=new byte[80] ;
        int i,j,k,l ;
        int ii,jj,kk ;
        int il,ih,jl,jh ;
        double alpha,beta,ff ;
        ppmread objImage = new ppmread() ;
        newwidth=960;
        newheight=720;
        if(args.length == 2) {
            objImage.ppmmain(args[0]) ;
            newwidth=Integer.parseInt(args[1]) ;
            newheight=(objImage.height*newwidth)/objImage.width ; 
        } else {
            System.err.println
                ("Usage: java ppmScaling filename new-width") ;
            System.exit(0) ;
        }
        // enlarge =1 for enlargement else enlarge=0
        if(newwidth>=objImage.width) enlargex=1 ;
        else enlargex=0 ;
        if(newheight>=objImage.height) enlargey=1 ;
        else enlargey=0 ;
        double xcoord[]=new double[objImage.width] ;
        double ycoord[]=new double[objImage.height] ;
        double newxcoord[]=new double[newwidth] ;
        double newycoord[]=new double[newheight] ;
        int imagData[][]=new int[newwidth][newheight] ;
        imageData=new short [newwidth][newheight] ;
        System.out.println("P"+objImage.magicidno) ;
        System.out.println("# Data is created by ppmImagSmooth") ;
        System.out.println("# filename = "+ objImage.sfilename ) ;
        System.out.println("# Width    = "+objImage.width+"-->"+newwidth) ;
        System.out.println("# Height   = "+objImage.height+"-->"+newheight) ;
        for(i=0;i<objImage.width;i++) xcoord[i]=(double)i/objImage.width ;
        for(i=0;i<objImage.height;i++) ycoord[i]=(double)i/objImage.height ;
        for(i=0;i<newwidth;i++) newxcoord[i]=(double)i/newwidth ;
        for(i=0;i<newheight;i++) newycoord[i]=(double)i/newheight ;
        //
        // Assign old coordinates for new coordinate calculation
        //
        int maxImage=0;
        for(i=0;i<newwidth;i++) {
            if (i==0) {
                alpha=0.0 ;
                il=0 ;
                ih=il+1 ;
            } else {
                if((enlargex==1) && (i == (newwidth-1))) {
                    ih=objImage.width-1 ;
                    il=ih;
                    alpha=0.0 ;
                } else {
                    for(ii=0;ii<objImage.width;ii++) {
                        if(xcoord[ii]>newxcoord[i]) break ;
                    } 
                    il =ii-1 ;
                    ih= ii ;
                    alpha=(newxcoord[i]-xcoord[il])/(xcoord[ih]-xcoord[il]) ;
                }
            }
            for(j=0;j<newheight;j++) {
                if (j==0) {
                    beta=0.0 ;
                    jl=0;
                    jh=jl+1;
                } else {
                    if((enlargey==1) && (j==(newheight-1))) {
                        beta=0.0 ;
                        jh=objImage.height-1 ;
                        jl=jh ;
                    } else {
                        for(jj=0;jj<objImage.width;jj++) {               
                            if(ycoord[jj]>newycoord[j]) break ;
                        }
                        jl=jj-1 ;
                        jh=jj ;
                        beta=(newycoord[j]-ycoord[jl])/(ycoord[jh]-ycoord[jl]) ;
                    }
                }
                ff =  (1.0-alpha)*(1.0-beta)*objImage.gdata[il][jl] ;
                ff += alpha*(1.0-beta)*objImage.gdata[ih][jl] ;
                ff += (1.0-alpha)*beta*objImage.gdata[il][jh] ;
                ff += alpha*beta*objImage.gdata[ih][jh] ;
                //System.out.println("("+i+","+j+") = "+ff) ;
                imageData[i][j]=(short)Math.round(ff) ;
                if(imageData[i][j] > maxImage) maxImage=imageData[i][j] ;
            }
        }

        if(objImage.magicidno>1) 
            System.out.println("# Maxval   = "+maxImage) ;
        System.out.println(newwidth+" "+newheight+" "+maxImage) ;

							
        int ndataMaxLine=50 ;
        int rdataline,rdatalast ;
        rdataline = newwidth/ndataMaxLine ;
        rdatalast = newwidth%ndataMaxLine ;
        for(k=0;k<newheight;k++) {
            for(i=0;i<rdataline;i++) {
                for(j=0;j<ndataMaxLine;j++) {
                    int ij=i*ndataMaxLine+j ;
                    System.out.print(imageData[ij][k]+" ") ;
                }
                System.out.println() ;
            }
            for(j=0;j<rdatalast;j++) {
                int ij=rdataline*ndataMaxLine+j ;
                System.out.print(imageData[ij][k]+" ") ;
            }
            System.out.println() ;
        }
    }
}

実行は

> java ppmScaling tmuinthouse2.pgm 720 > inthouse2720.pgm

などとする.オリジナルの縦横比は保存されるようにプログラムで設定している.オリジナルのピクセル幅×高さは480×360なので1.5倍の拡大となる.

tmuinthouse2.pgm

inthouse2720.pgm

移動・回転

アフィン変換

課題3

拡大縮小・移動・回転・反転を以下のひとつの式で行うアフィン変換のプログラムを作成しそのどれにも対応できるよう(例えば拡大縮小を1、移動を2、回転を3、反転を4、5としてさらに拡大縮小+移動あるいは拡大縮小+回転あるいは拡大縮小+移動+回転さらに反転も加えてそれらにそれぞれ数値を与えて)プログラムの初期パラメータを設定し原図ファイルを思うままに変換させるプログラムを作成せよ.

アフィン変換は以下のように与えられる;

それぞれの係数の意味は以下のように与えられる.

なお上記の式、図は愛媛大学、村上・泉田研究室の「画像処理」より引用しました.

2値化処理

階層均一化



C.Fukunaga