Julia でヒストグラム ― 2016年02月02日 19:30
話は前後するが、Julia でヒストグラム表示をやってみたのでメモ
イメージの読み出しが load() になった。読み出すと Image 構造体に読み込まれるので、その中から data() で画像データのみを取り出す。取り出した data は、RGB 構造体の配列なので、reinterpret () で UInt8 の配列に変換する。これで、3 x 横 x 縦 サイズの UInt8 の配列に取り出せる。
さらに、RGB に分離する。前記 UInt8 配列のうち、1次元目はR/G/B をあらわしているので、R のみ取り出すには [1,:,:] を指定する。このままでは3次元の配列なので2次元にするために、先頭の次元を squeze で削除する。これでやっと2次元の R が取り出せる。(結構面倒くさい)
続いてヒストグラム生成。2次元のままだと(所望の)ヒストグラムが取れないので、reshape で1次元に並び替える。squeeze とは違うので注意(前述の次元削除も reshape で出来ないこともないけど)。
これでやっとヒストグラムが取れる。hist(v, e) の書式で、v はヒストグラムを取るベクトル(1次元の配列)。e は幅指定で、-1:255 で 0から255の範囲をヒストグラム化。1番目の範囲は、-1より大きく、0以下。255番目の範囲は245より大きく、255以下が含まれる。
ちなみに、-1:16:255 でヒストグラムの幅を16ずつとかに出来る。edge と 分布が返ってくる。
ヒストグラムが出来たので、グラフに描いてみる。ここでは、描画パッケージとして Winston を使用。
Histogram() でグラフのオブジェクトを作成、FramedPlot() は、描画用のウィンドウの作成。add() で、ウィンドウにグラフオブジェクトを追加。これでヒストグラムが表示される。
やってみてわかったのだが、UInt8 のまま hist() を呼ぶと動作がおかしい。* 1.0 で Float64 に変換しているのはそのため。理由は不明・・・
UInt8 のままで画像処理をするといろいろ面倒なので、最初に RGB 変換時に Float に変換したほうがいいかも知れない。
イメージの読み出しが load() になった。読み出すと Image 構造体に読み込まれるので、その中から data() で画像データのみを取り出す。取り出した data は、RGB 構造体の配列なので、reinterpret () で UInt8 の配列に変換する。これで、3 x 横 x 縦 サイズの UInt8 の配列に取り出せる。
さらに、RGB に分離する。前記 UInt8 配列のうち、1次元目はR/G/B をあらわしているので、R のみ取り出すには [1,:,:] を指定する。このままでは3次元の配列なので2次元にするために、先頭の次元を squeze で削除する。これでやっと2次元の R が取り出せる。(結構面倒くさい)
続いてヒストグラム生成。2次元のままだと(所望の)ヒストグラムが取れないので、reshape で1次元に並び替える。squeeze とは違うので注意(前述の次元削除も reshape で出来ないこともないけど)。
これでやっとヒストグラムが取れる。hist(v, e) の書式で、v はヒストグラムを取るベクトル(1次元の配列)。e は幅指定で、-1:255 で 0から255の範囲をヒストグラム化。1番目の範囲は、-1より大きく、0以下。255番目の範囲は245より大きく、255以下が含まれる。
ちなみに、-1:16:255 でヒストグラムの幅を16ずつとかに出来る。edge と 分布が返ってくる。
ヒストグラムが出来たので、グラフに描いてみる。ここでは、描画パッケージとして Winston を使用。
Histogram() でグラフのオブジェクトを作成、FramedPlot() は、描画用のウィンドウの作成。add() で、ウィンドウにグラフオブジェクトを追加。これでヒストグラムが表示される。
やってみてわかったのだが、UInt8 のまま hist() を呼ぶと動作がおかしい。* 1.0 で Float64 に変換しているのはそのため。理由は不明・・・
UInt8 のままで画像処理をするといろいろ面倒なので、最初に RGB 変換時に Float に変換したほうがいいかも知れない。
# image 読み込み
using Images
img= reinterpret(UInt8,data(load("SleepingCat.bmp")))
red = squeeze(img[1,:,:], 1)
green= squeeze(img[2,:,:], 1)
blue = squeeze(img[3,:,:], 1)
# Histogram 生成
data1d= reshape(green, length(green)) * 1.0
(edges, counts)= hist(data1d, -1:255)
# Winston の機能
using Winston
h= Histogram(edges, counts)
p= FramedPlot()
add(p, h)
julia - bicubic ― 2015年06月18日 19:51
バイキュービックの計算には、16画素のデータそれぞれに対する係数を求める必要がある。
求める画素位置と、実際の画素の位置から、画素間の距離を求め、その距離に応じた係数をそれぞれの画素の輝度値に掛け、その総和を求める画素の輝度値とする。
ここで、距離と係数の関係は sinc 関数を3次近似したものが使われる。
上の画像のピンクの線が sinc 関数で、下の緑と青が近似した値である。見やすくするために sinc 関数は上に0.5だけ移動してある。
ちなみに、sinc 関数とは sin(πx)/(πx) のこと。
この部分を関数にするとこんな感じ。
x方向、y方向、それぞれ独立に係数を求め、掛け合わせることで16個の係数を求める。(試しに xy両方を含めた距離で計算してみたが、ノイズっぽい画像になってしまった)
x方向だけ考えると、求める画素の位置は、4つの画素の真ん中にくる。
x値を小さい順に x1, x2, x3, x4 と置くと、求める画素の位置(x)は x2≦x<x3 となり、x2= floor(x) である。
よって、x と x2 の距離は x - floor(x) であり、x と x1 の距離は x と x2 の距離に1を足したものである。
また、x と x3 の距離は、1から(x と x2 の距離)を引いたものとなり、x と x4 の距離はそれに1を足したものである。
よって、x と x2 の距離(すなわち、x - floor(x) )を xd とすれば、それぞれ xd+1, xd, 1-xd, 2-xd となる。
y に関しても同様。
これを関数で書くと以下のようになる。
map を使って、行列の距離値を sinc 関数を使った係数に変換している。
xc は、4行1列の行列。yc は、1行4列の行列。
これらを掛け合わせると、4行4列の行列が得られる。これを戻り値にしている。
完成した関数は以下のとおり。
求める画素位置と、実際の画素の位置から、画素間の距離を求め、その距離に応じた係数をそれぞれの画素の輝度値に掛け、その総和を求める画素の輝度値とする。
ここで、距離と係数の関係は sinc 関数を3次近似したものが使われる。
上の画像のピンクの線が sinc 関数で、下の緑と青が近似した値である。見やすくするために sinc 関数は上に0.5だけ移動してある。
ちなみに、sinc 関数とは sin(πx)/(πx) のこと。
この部分を関数にするとこんな感じ。
function fsinc(d)
return (d < 1) ? d^3 - 2d^2 + 1 : (d < 2) ? -d^3 + 5d^2 - 8d + 4 : 0
end
実は julia には sinc 関数があるので、別の名前にした。x方向、y方向、それぞれ独立に係数を求め、掛け合わせることで16個の係数を求める。(試しに xy両方を含めた距離で計算してみたが、ノイズっぽい画像になってしまった)
x方向だけ考えると、求める画素の位置は、4つの画素の真ん中にくる。
x値を小さい順に x1, x2, x3, x4 と置くと、求める画素の位置(x)は x2≦x<x3 となり、x2= floor(x) である。
よって、x と x2 の距離は x - floor(x) であり、x と x1 の距離は x と x2 の距離に1を足したものである。
また、x と x3 の距離は、1から(x と x2 の距離)を引いたものとなり、x と x4 の距離はそれに1を足したものである。
よって、x と x2 の距離(すなわち、x - floor(x) )を xd とすれば、それぞれ xd+1, xd, 1-xd, 2-xd となる。
y に関しても同様。
これを関数で書くと以下のようになる。
function getCOEFFS (xin, yin)
# coefficient (a= -1)
function fsinc(d)
return (d < 1) ? d^3 - 2d^2 + 1 : (d < 2) ? -d^3 + 5d^2 - 8d + 4 : 0
end
# calc coefficients
xd= xin - floor(xin)
xc= map(x -> fsinc(x), [xd+1; xd; 1-xd; 2-xd])
#
yd= yin - floor(yin)
yc= map(x -> fsinc(x), [yd+1 yd 1-yd 2-yd])
#
return xc * yc
end
関数の中に関数が書けるのは便利。map を使って、行列の距離値を sinc 関数を使った係数に変換している。
xc は、4行1列の行列。yc は、1行4列の行列。
これらを掛け合わせると、4行4列の行列が得られる。これを戻り値にしている。
完成した関数は以下のとおり。
# ibuff の (xin,yin) から画素値を取得する
function imgpic (ibuff, xin, yin)
# get 4x4 pixels
function getPIX4x4 (ibuff, xin, yin)
xstart= ifloor(xin) - 1
ystart= ifloor(yin) - 1
try
return ibuff[xstart:(xstart+3), ystart:(ystart+3)]
catch # out of range
obuff= ibuff[1:4,1:4]
(xmax,ymax)= size(ibuff)
for x= 1:4, y= 1:4
xin= xstart + (x - 1)
yin= ystart + (y - 1)
xin= (xin < 1) ? 1 : (xin > xmax) ? xmax : xin
yin= (yin < 1) ? 1 : (yin > ymax) ? ymax : yin
obuff[x,y]= ibuff[xin,yin]
end
return obuff
end
end
# get 4x4 coefficients
function getCOEFFS (xin, yin)
# coefficient (a= -1)
function fsinc(d)
return (d < 1) ? d^3 - 2d^2 + 1 : (d < 2) ? -d^3 + 5d^2 - 8d + 4 : 0
end
# calc coefficients
xd= xin - floor(xin)
xc= map(x -> fsinc(x), [xd+1; xd; 1-xd; 2-xd])
#
yd= yin - floor(yin)
yc= map(x -> fsinc(x), [yd+1 yd 1-yd 2-yd])
#
return xc * yc
end
pixels= getPIX4x4(ibuff, xin, yin)
coeffs= getCOEFFS(xin, yin)
return sum(pixels .* coeffs) / sum(coeffs)
end
julia - 拡大・縮小 ― 2015年06月12日 18:14
フィルター処理だけであれば、imfilter で簡単に行うことが出来る。
今回は、画像の拡大・縮小の処理を行ってみる。
最終的には、imgzoom( img, xsize, ysize) という関数で、画像を拡大・縮小できるようにしたい。
任意サイズに画像を拡大・縮小する場合、画素の補間が必要となる。
ここでは、bicubic 法で補完を行う。
バイキュービックでは、周辺の16画素を用いて補完を行う。この際には、補完係数に3次関数を用いるので cubic と呼ばれるようだ。ちなみに、この補完係数は、sinc関数を3次近似したものらしい。
ということで、(x,y) の周辺16画素を取得する関数をまずは作ってみる。
簡単に書けばこんな感じか。
しかし、これを毎回行うのは効率的ではないので、エラーが発生した時だけ真面目に処理をするようにした。(ちょっとだけ早い)
とりあえず今日はここまで。
今回は、画像の拡大・縮小の処理を行ってみる。
最終的には、imgzoom( img, xsize, ysize) という関数で、画像を拡大・縮小できるようにしたい。
任意サイズに画像を拡大・縮小する場合、画素の補間が必要となる。
ここでは、bicubic 法で補完を行う。
バイキュービックでは、周辺の16画素を用いて補完を行う。この際には、補完係数に3次関数を用いるので cubic と呼ばれるようだ。ちなみに、この補完係数は、sinc関数を3次近似したものらしい。
ということで、(x,y) の周辺16画素を取得する関数をまずは作ってみる。
簡単に書けばこんな感じか。
function getPIX4x4 (img, x, y)
xstart= ifloor(x) - 1
ystart= ifloor(y) - 1
return img[xstart:(xstart+3), ystart:(ystart+3)]
end
ただし、これだと、x=1 の時(左端)や右端の時に img の添え字をはみ出してしまうので、x, y の値のチェックが必要。しかし、これを毎回行うのは効率的ではないので、エラーが発生した時だけ真面目に処理をするようにした。(ちょっとだけ早い)
function getPIX4x4 (img, xin, yin)
xstart= ifloor(xin) - 1
ystart= ifloor(yin) - 1
try
return img[xstart:(xstart+3), ystart:(ystart+3)]
catch
buff= img[1:4,1:4]
(xmax,ymax)= size(img)
for x= 1:4, y= 1:4
xpos= xstart + (x - 1)
ypos= ystart + (y - 1)
xpos= (xpos < 1) ? 1 : (xpos > xmax) ? xmax : xpos
ypos= (ypos < 1) ? 1 : (ypos > ymax) ? ymax : ypos
buff[x,y]= img[xpos,ypos]
end
return buff
end
end
try ブロックを処理中にエラーが生じると catch ブロックに処理が移る。とりあえず今日はここまで。
julia - BMP の読み書き ― 2015年06月05日 17:36
Images のパッケージが大幅に変更されました。
imread(), imwrite() は、それぞれ load(), save() に変わりました。
機会があれば対応したものを執筆しますが、混乱を避けるため、とりあえず消しておきます。(2015/12/19)
画像処理をするために、BMPファイルの読み書きをやってみた。
Pkg.add("Images") でパッケージを追加。
入出力には ImageMagick を使用しているらしい。
以下、ちょっとしたコード。
RGB Image 構造体
ちなみに、data は省略できて、img1[100,100].r でも同じ。
ただし、img1[:,:].r のように範囲指定はできない。
構造体の properties は、img1["IMcs"] のように指定する。
RGB それぞれを行列として扱いたい場合には、red(), green(), blue() という関数で分離できる
imread で読み込んだ場合、この行列は、各画素 0 から 1 の値となる。
255 を掛けると元の BMP のピクセル値に一致する。
Image 構造体を作る
行列から Image 構造体を作るには、grayim(),colorim()関数を使う。
grayim() は、m x n の行列からグレースケールの Gray Image 構造体を作る。
colorim() は、3 x m x n の行列が必要。RGB から直接変換する方法は不明。
とりあえず、この大きさの行列を作って、RGB 値を代入すれば使える。
まとめ
imread(), imwrite() は、それぞれ load(), save() に変わりました。
機会があれば対応したものを執筆しますが、混乱を避けるため、とりあえず消しておきます。(2015/12/19)
Pkg.add("Images") でパッケージを追加。
入出力には ImageMagick を使用しているらしい。
以下、ちょっとしたコード。
# Images パッケージを使用する宣言
using Images
# sample.bmp を読み込む
img1= imread("sample.bmp")
# 白黒画像に変換
img1= convert(Image{Gray}, img1)
# sobel filter の係数を取得。2つの行列を返してくる
sobelX, sobelY= sobel()
# X方向の sobel filter を適用
img1= imfilter(img1, sobelX)
# 結果をファイルに書き込む
imwrite(colorim(img1), "edgeX.bmp")
RGB Image 構造体
img1= imread("sample.bmp")
評価値
RGB Image with:
data: 1280x720 Array{RGB{UfixedBase{Uint8,8}},2}
properties:
IMcs: sRGB
spatialorder: x y
pixelspacing: 1 1
ピクセルデータには、img1.data[100,100].r のようにしてアクセスできる。ちなみに、data は省略できて、img1[100,100].r でも同じ。
ただし、img1[:,:].r のように範囲指定はできない。
構造体の properties は、img1["IMcs"] のように指定する。
RGB それぞれを行列として扱いたい場合には、red(), green(), blue() という関数で分離できる
imgR= red (img1)
imgG= green(img1)
imgB= blue (img1)
これで imgR,imgG,imgB は、行列となるので、普通に計算などができる。imread で読み込んだ場合、この行列は、各画素 0 から 1 の値となる。
255 を掛けると元の BMP のピクセル値に一致する。
Image 構造体を作る
行列から Image 構造体を作るには、grayim(),colorim()関数を使う。
grayim() は、m x n の行列からグレースケールの Gray Image 構造体を作る。
colorim() は、3 x m x n の行列が必要。RGB から直接変換する方法は不明。
とりあえず、この大きさの行列を作って、RGB 値を代入すれば使える。
まとめ
using Images
# sample.bmp 読み出し
img1= imread("sample.bmp")
# RGB に分離
imgR= red (img1)
imgG= green(img1)
imgB= blue (img1)
# Gray を作ってみる
imgW= imgR * 0.29 + imgG * 0.60 + imgB * 0.11
imwrite(grayim(imgW), "gray.bmp")
# RGB Image 構造体に戻す
imgRGB= zeros(3, size(imgR,1), size(imgR,2))
imgRGB[1,:,:]= imgR
imgRGB[2,:,:]= imgG
imgRGB[3,:,:]= imgB
imwrite(colorim(imgRGB), "rgb.bmp")
最近のコメント