（す）メモ帳

バイキュービックの計算には、16画素のデータそれぞれに対する係数を求める必要がある。

求める画素位置と、実際の画素の位置から、画素間の距離を求め、その距離に応じた係数をそれぞれの画素の輝度値に掛け、その総和を求める画素の輝度値とする。
ここで、距離と係数の関係は sinc 関数を３次近似したものが使われる。
上の画像のピンクの線が 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方向だけ考えると、求める画素の位置は、４つの画素の真ん中にくる。
x値を小さい順に x1, x2, x3, x4 と置くと、求める画素の位置(x)は x2≦x＜x3 となり、x2= floor(x) である。
よって、x と x2 の距離は x - floor(x) であり、x と x1 の距離は x と x2 の距離に１を足したものである。
また、x と x3 の距離は、１から(x と x2 の距離)を引いたものとなり、x と x4 の距離はそれに１を足したものである。
よって、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行１列の行列。yc は、１行４列の行列。
これらを掛け合わせると、４行４列の行列が得られる。これを戻り値にしている。

完成した関数は以下のとおり。

# 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