画像解析(演習6)

by 菊地時夫 — last modified 2008-11-10 08:21

~/image/e6 で演習

ひとくちに

  • 画像解析と言ってもピンからキリまで。
  • 何を明らかにしたいのかによって方法は様々。
  • なので、今回は「画像解析ごっこ」

とりあえず、クラスタリング

  • だが、階層的クラスタリングは計算時間がかかるので、今日はやらない
  • 非階層的クラスタリングを使う
  • 計算を走らせるたびに違った結果になったりもする。

サンプル画像

  • 以前使用した 1.ppm を使う。
  • ダウンロードして本日の演習ディレクトリに入れておく

共通で使用する後処理

  • クラスタリングした結果を見るのに、
  • そのクラスタに属する点の平均(R,G,Bごと)を計算して
  • 画像で表示(または保存)する

ClusterImg.py

# -*- coding: utf-8 -*-
# ClusterImg.py
import Image
 
def calcavr(clusterid, data, N):
    """ クラスタごとに (R,G,B) の平均を求める
    clusterid ... クラスタリング結果
    data ........ 元のデータ (R,G,B) のリスト
    N ........... クラスタ数
"""
    clustavr = []
    for i in range(N):
        R = G = B = 0
        count = 0
        pos = 0
        for x in list(clusterid):
            if x == i:
                r, g, b = data[pos]
                R += r
                G += g
                B += b
                count += 1
            pos += 1
        R = R/count
        G = G/count
        B = B/count
        clustavr.append((R,G,B))
    return clustavr
 
def clusterimg(size, clusterid, clustavr):
    """ クラスタリングの結果をそれぞれのクラスタ内の平均値で色を付けて画像にする
    size ... 画像サイズ (width, height) のタプル
    clusterid ... クラスタリング結果
    clustavr .... 各クラスタの平均値(R,G,B)
"""
    cimg = Image.new('RGB', size)
    cdata = [clustavr[clusterid[i]] for i in range(size[0]*size[1])]
    cimg.putdata(cdata)
    return cimg

そのままクラスタリング

  • できたらプログラムを修正して、結果の画像を s601.ppm に保存する。
    # -*- coding: utf-8 -*-
    # e601.py
    import Image
    from Pycluster import *
    import ClusterImg
     
    img = Image.open('1.ppm')
    width , height = img.size
    N = 10
    data = list(img.getdata())
    # kcluster をデフォルト設定で使用
    clusterid, error, nfound = kcluster(data, N)
    # 後処理と表示
    clustavr = ClusterImg.calcavr(clusterid, data, N)
    print clustavr
    cimg = ClusterImg.clusterimg(img.size, clusterid, clustavr)
    cimg.show()
    

前処理で情報を追加

  • 次のような処理をクラスタリングの前に入れると、高波数成分が多いという情報を 追加できる(はず)
    imx = img.convert('L')
    imx = imx.filter(ImageFilter.FIND_EDGES)
    imx = imx.filter(ImageFilter.SMOOTH_MORE)
    enh = ImageEnhance.Brightness(imx)
    imx = enh.enhance(6.0)
    datx = Numeric.transpose(Numeric.array([list(imx.getdata())]))
    data = Numeric.array(img.getdata())
    datx = Numeric.concatenate((data, datx), 1)
    
  • ImageFilter, ImageEnhance, Numeric など必要なモジュールを import しておくことも忘れずに。
  • プログラムは e602.py で作成、結果は s602.ppm に保存。

前処理で情報を減らす

  • 画像の濃淡情報を削除して、色情報だけでクラスタリングしてみる
  • YUV 変換で UV 成分だけを使用する。
    uv = []
    for x in data:
        r, g, b = x
        u = -0.14713 * r - 0.28886 * g + 0.436   * b
        v =  0.615   * r - 0.51499 * g - 0.10001 * b
        uv.append((u, v))
    
  • この uv をクラスタリングするための data として与える。
  • プログラムを e603.py で作成、結果は s603.ppm に保存。

突然ですが、主成分分析

  • Principal Component Analysis
  • 共分散行列の固有値問題の解
  • 主成分分析によって得られた成分でデータを表す
  • =特異値分解(と、マニュアル (p.29) に書いてある)

主成分分析の練習

  • Gnuplot を使うので X11 で実行
    # e604.py ... PCA example
    import Numeric
    import LinearAlgebra
    import random
    import Gnuplot
     
    # Generate skewed random data.
    z = []
    for i in range(100):
        x = random.gauss(0, 2.0)
        y = random.gauss(0, 1.0)
        z.append(((x + y), (x - y)))
     
    z = Numeric.array(z)
    g1 = Gnuplot.Gnuplot()
    g1.plot(z)
    # Perform PCA
    a = LinearAlgebra.singular_value_decomposition(z)
    g2 = Gnuplot.Gnuplot()
    g2.plot(a[0])
    raw_input('done? ')
    
  • 2つのプロットデータがどのように分布しているか違いを見る

発展課題ではあるが解答付き

  • 気象情報頁 で公開している画像から 2007年7月6日03UTC の画像を選び、 日本付近を切り出しました。
  • IR1.pgm , IR2.pgm , IR3.pgm , IR4.pgm , VIS.pgm (ダウンロードして、本日の作業ディレクトリに入れてください)
  • 主成分分析の第3成分までを使ってクラスタリングする。
  • ただし、色分けはクラスタごとに白、黒、赤、緑、青、イエロー、シアン、マゼンタを使う。

プログラム

# e605.py
import Image
import Numeric
import LinearAlgebra
from Pycluster import *
 
img1 = Image.open('IR1.pgm')
img2 = Image.open('IR2.pgm')
img3 = Image.open('IR3.pgm')
img4 = Image.open('IR4.pgm')
img5 = Image.open('VIS.pgm')
width , height = img1.size
N = 8
 
dat1 = Numeric.transpose(Numeric.array([list(img1.getdata())]))
dat2 = Numeric.transpose(Numeric.array([list(img2.getdata())]))
dat3 = Numeric.transpose(Numeric.array([list(img3.getdata())]))
dat4 = Numeric.transpose(Numeric.array([list(img4.getdata())]))
dat5 = Numeric.transpose(Numeric.array([list(img5.getdata())]))
data = Numeric.concatenate((dat1, dat2, dat3, dat4, dat5), 1)
x = LinearAlgebra.singular_value_decomposition(data)
y = [(p[0], p[1], p[2]) for p in x[0]]
clusterid, error, nfound = kcluster(y, N)
 
clustclr = [(0,0,0),(255,0,0),(0,255,0),(0,0,255),
            (255,255,0),(0,255,255),(255,0,255),(255,255,255)]
cimg = Image.new('RGB', img1.size)
cdata = [clustclr[clusterid[i]] for i in range(len(data))]
cimg.putdata(cdata)
cimg.show()

できたら

  • s605.gif に画像を保存。(クラスタリングの結果は毎回違うと思うが)

まいど、出席確認

  • 本日作成した 最後のプログラム をメールの本文に記入して、菊地・画像処理論用 に送信。本文の最初に感想等入れるも可。