카테고리 없음
[파이썬/Ctypes] Ctypes을 사용해 4차원 합성곱 구하는 함수 (최종)
미친토끼
2021. 3. 12. 12:17
2차원 합성곱을 C함수로 구하는 Ctypes 기능을 4차원 파이썬 코드에 붙인 것이다.
이틀 전에 만든 directConvolution.py에서 제일 말단 함수인 directConvolution3D1F의 반복문을 수정해서
directConvolve2d.c를 사용할 수 있게 했다.
np.zeros 사용하는 부분에 영어로 길게 주석을 붙여놨는데, 결론은, 변수를 재사용할 가능성이 있을 때에는 해당 변수를 초기화(여기서는 0으로)하는 게 정신건강에 좋다는 것이다. 이 코드를 테스트하면서 out을 재사용하고 있는데 이상하게 합계가 누적되는 경향이 있어서 조사해보니 np.zeros로 한번 0으로 초기화해둔 것을 믿고, 두 번 다시 초기화하지 않았다는 데 있었다. 그래서 C코드에서 아예 0으로 초기화를 했다.
# This code was written by Don Han, 2021.3.12
# https://blog.naver.com/madrabbit7
# filename: directConvolutionCtypes.py
# Sorry for my poor english comments.^^
# My (f***ing) spyder IDE in fedora 33 didn't allow me to write Haugul here.
# some of hangul comments came from another IDE (pyCharm or Studio Code etc...)
# This is a direct convolution function which impliments the definition
# of original convolution concept.
# First, I referred a code of some 3 dimenstional direct convolution function here:
# https://welcome-to-dewy-world.tistory.com/94
# So I made this 4 dimensional version which has modulized functions.
# Special thanks to Mr. Han Kyeonghoon, the professor of Data Science department of Suwon Univ.
# He gave me a lot of knowledges and inspirations.
""" Function List
pad_with(vector, pad_width, iaxis, kwargs):
numpy pad function
directConvolution3D1FCtypes(x, f, b=0, stride=1):
3D convolution function, accept only 1 filter, accept 1 data
def directConvolution3DNF(x, f, b=0, stride=1):
3D convolution function. accept N filters, accept 1 data.
It calls directConvolution3D1F function.
directConvolution(x, f, b=0, stride=1, pad=0):
4D convolution function. accept N data, N filters.
It does padding job, checks bias array, and has a main loop.
It calls directConvolution3DNF function.
"""
import numpy as np
import ctypes
import sys
import os.path
np.set_printoptions(precision=0, suppress=True)
# Numpy pad Fuction
# I got this code from https://numpy.org/doc/stable/reference/generated/numpy.pad.html
def pad_with(vector, pad_width, iaxis, kwargs):
"""
Usage: padded_x = np.pad(x, 2, pad_with)
Be careful! np.pad transformed the shape of x to another.
So you don't broadcast like this:
x = np.pad(x, 2, pad_with)
if so, immediately your computer will be exploded, next moments,
you will find yourself be in Mars with your precious Perseverance forever!
"""
pad_value = kwargs.get('padder', 0)
vector[:pad_width[0]] = pad_value
vector[-pad_width[1]:] = pad_value
def directConvolution3D1FCtypes(x, f, b=0, stride=1):
"""
It implements original convolution concept. It can be operated independently.
x: input data (3 dimension only)
f: filter (3 dimension only, and 1 filter only)
b: bias
"""
if x.ndim != 3 or f.ndim != 3:
sys.exit(f"The data({x.ndim}) and filter dimensions({f.ndim}) all should be 3.")
C, xh, xw = x.shape
Cf, fh, fw = f.shape
if C != Cf:
sys.exit(f"The number of channel of data({C}) and the number of \
channel of filters({Cf}) should be same.")
x = np.float64(x)
f = np.float64(f)
oh = 1 + int((xh-fh)/stride)
ow = 1 + int((xw-fw)/stride)
# In the C function, 'out' arrays again will be initialized to 0.
# because in python, we can reuse 'out' arrays by mistake
# without initializing to 0.
# Firstly I believed below np.zeros's initializing, and reused 'out' arrays
# without initialing to 0, then whenever Ctypes got return value of 'out' from
# the C function, I found out that in 'out' arrays, their values had been
# acculmulated by then. So I added 'initializing to 0' code to the C function like this,
# dst[i * ow + j] = 0;
# Of course, you can initialize the value in Python code also before reusing them.
out = np.zeros((oh, ow), dtype=np.float64)
# This valiable is for sum of each channels.
total = out.copy()
c_lib = ctypes.CDLL("./directConvolve2d.so")
for i in range(C): # 채널 개수만큼 반복
src = x[i].ctypes.data_as(ctypes.POINTER(ctypes.c_double))
filter = f[i].ctypes.data_as(ctypes.POINTER(ctypes.c_double))
# I reuse 'out' arrays (C) times. C function code will initialize the value.
dest = out.ctypes.data_as(ctypes.POINTER(ctypes.c_double))
c_lib._directConvolve2d(src, xh, xw, filter, fh, fw, dest, stride)
total += out
return total + b
def directConvolution3DNF(x, f, b=0, stride=1):
"""
It implements original convolution concept. It can be operated independently.
It does not have padding module. Main function has that.
x: input data (3 dimension only)
f: filter (4 dimension available, N filters)
b: bias
"""
if x.ndim != 3 or f.ndim != 4:
sys.exit("Input data dimension sould be 3, and filter dimension should be 4.")
xc, xh, xw = x.shape
fn, fc, fh, fw = f.shape
if xc != fc:
sys.exit(f"Number of data channels({xc}) and filter channels({fc}) should be same!")
oh = int((xh - fh)/stride + 1)
ow = int((xw - fw)/stride + 1)
output = np.zeros([oh, ow])
# if b is an integer or a float number, then we make it an arrary of the same number
if type(b) == int or type(b) == float:
b = np.full((fn,), b)
for i in range(fn):
out = directConvolution3D1FCtypes(x, f[i], b[i], stride)
if i == 0:
link = out
else :
link = np.vstack([link, out])
return link.reshape(fn, oh, ow)
def directConvolution(x, f, b=0, stride=1, pad=0):
"""
x: input data (4 dimension only)
f: filter (4 dimension only, N filters)
b: bias
"""
xn, xc, xh, xw = x.shape
fn, fc, fh, fw = f.shape
if xc != fc:
sys.exit("A number of input data's channels should be same as the number of filter channels.")
if x.ndim != 4 or f.ndim != 4:
string = f"The dimensions of input data({x.ndim}) and filters({f.ndim}) should be 4 \
dimensions."
sys.exit(string)
# padding job
if type(pad) == int and pad > 0:
# I made new array variable which got the padded array.
padded_x = np.zeros((xn, xc, xh + 2*pad, xw + 2*pad))
for i in range(xn): # repeat as number of input data
for j in range(xc): # repeat as number of input data's chennel.
# The padded matrix can't be broadcasted to the same matrix.
# because of different height and width. for example,
# x[i][j] = np.pad(x[i][j], pad, pad_with)
# wil issue some error about broadcast.
padded_x[i][j] = np.pad(x[i][j], pad, pad_with)
elif pad == 0:
padded_x = x
# again, because, this is after padding
xn, xc, xh, xw = x.shape
fn, fc, fh, fw = f.shape
oh = int((xh - fh)/stride + 1)
ow = int((xw - fw)/stride + 1)
# if b is an integer or a float number, then we make it to an arrary of the same number.
if type(b) == int or type(b) == float:
b = np.full((fn,), b)
elif type(b) == str:
sys.exit("You shouldn't give a string as the bias.")
elif len(b) != fn:
sys.exit("A number of bias should be same as the number of filters.")
# this is the main loop
for i in range(xn): # xn = the number of input data
out = directConvolution3DNF(padded_x[i], f, b, stride)
if i == 0:
link = out
else :
link = np.vstack([link, out])
#link = np.hstack([link, out]) :same results ?
return link
if __name__ == "__main__":
x = np.array([[[1, 0, 1], [0, 1, 0], [1, 0, 1]],
[[1, 1, 1], [1, 0, 1], [1, 1, 1]]])
x2 = np.array([[[[1, 1, 1], [1, 0, 0], [1, 0, 0]],
[[1, 1, 1], [0, 0, 1], [0, 0, 1]]],
[[[1, 0, 0], [1, 0, 0], [1, 1, 1]],
[[0, 0, 1], [0, 0, 1], [1, 1, 1]]]])
f = np.array([[[1, 0], [1, 0]],
[[1, 1], [0, 0]]])
f2 = np.array([[[[1, 0], [1, 0]],
[[1, 1], [0, 0]]],
[[[0, 1], [0, 1]],
[[0, 0], [1, 1]]]])
# out = directConvolution3DNF(x, f2, b=[1,10])
# print("-----------------\n", out)
# (3, 4, 4)
x3=np.array([[[1,2,3,0],
[0,1,2,3], # --> red 행렬
[3,0,1,2],
[2,3,0,1]],
[[2,3,4,1],
[1,2,3,4], # --> green 행렬
[4,1,2,3],
[3,4,1,2]],
[[3,4,5,2], # --> blue 행렬
[2,3,4,5],
[5,2,3,4],
[4,5,2,3]]])
# (3, 3, 3)
f3=np.array([[[2,0,1],
[0,1,2],
[1,0,2]],
[[3,1,2],
[1,2,3],
[2,1,3]],
[[4,2,3],
[2,3,4],
[3,2,4]]])
#out = directConvolve3D1FCtypes(x3, f3, b=0, stride=1)
#print("-----------------\n", out)
# out = directConvolution(x2, f2, b=0, stride=2, pad=2)
# print(out)
# x4 = np.random.randint(0, 256, 256*256*30)
# x4 = x4.reshape((10, 3, 256, 256))
# f4 = np.random.randint(-2, 3, 81)
# f4 = f4.reshape(3, 3, 3, 3)
# b = [0, 2, 1]
# out = directConvolution(x4, f4, b=0, stride=1, pad=0)
# print(out[:100])
# We will use random datasets for test. For data verifying from this function \
# and the book's Convolution function, we need same datasets.
# So we use Numpy file save and load function.
# We make random datasets just once.
filename_data = "random_data.npy"
filename_filter = "random_filter.npy"
# if dataset files do not exist, then make them.
if os.path.isfile(filename_data) != True:
x4 = np.random.randint(0, 256, 100*3*256*256)
x4 = x4.reshape((100, 3, 256, 256))
np.save(filename_data, x4)
x4 = np.load(filename_data)
if os.path.isfile(filename_filter) == False:
f4 = np.random.randint(-2, 3, 81)
f4 = f4.reshape(3, 3, 3, 3)
np.save(filename_filter, f4)
f4 = np.load(filename_filter)
out = directConvolution(x4, f4, b=[1,2,0], stride=1, pad=0)
print(out[:2])