【iPhone イメージセンサ 解析】Sony製IMX503を徹底調査してみよう

import numpy as np
import matplotlib.pyplot as plt
import math

#***********************************
#1.画像設定
#***********************************
width=4032
height=3024
width_roi = 100
height_roi = 100
image_num = 15
bit = 12
dark = 528 #

half_width = int(width/2)
half_height = int(height/2)
half_width_roi = int(width_roi/2)
half_height_roi = int(height_roi/2)

#***********************************
#2.各種配列作成
#***********************************
img_roi=np.empty((height_roi, width_roi),dtype = np.uint16)
img_s = np.empty((4,image_num),dtype = np.float32)
img_n = np.empty((4,image_num),dtype = np.float32)
img_sn = np.empty((4,image_num),dtype = np.float32)
rgb = np.array(range(4)) 
rgb_det = rgb.reshape((2,2))
img_roi_r=np.empty((half_height_roi, half_width_roi),dtype = np.uint16)
img_roi_gr=np.empty((half_height_roi, half_width_roi),dtype = np.uint16)
img_roi_gb=np.empty((half_height_roi, half_width_roi),dtype = np.uint16)
img_roi_b=np.empty((half_height_roi, half_width_roi),dtype = np.uint16)
x_axis = np.empty((1,image_num),dtype = np.uint16)

#***********************************
#3.複数枚のRaw処理
#***********************************
for i in range(image_num):
    x_axis[0][i] = 2**i
    rawpath = "DNG_4032x3024_dark=528_12bit_lin"+str(i)+".raw"
    #***********************************
    #raw画像を開く
    #***********************************
    fd = open(rawpath, 'rb')
    f = np.fromfile(fd, dtype=np.uint16, count=height*width)
    img = f.reshape((height,width))
    fd.close()

    #***********************************
    #ROI切り出し
    #***********************************
    width_roi_start = int((width-width_roi)/2)
    height_roi_start = int((height-height_roi)/2)
    for x in range(width_roi_start,width_roi_start+width_roi):
        for y in range(height_roi_start,height_roi_start+height_roi):
            img_roi[y-height_roi_start][x-width_roi_start] = img[y][x]

    #***********************************
    #ダーク除去
    #***********************************
    for x in range(width_roi):
        for y in range(height_roi):
            if img_roi[y][x] > dark:
                img_roi[y][x] = img_roi[y][x]-dark
            else:
                img_roi[y][x] = 0

    #***********************************
    #SN計算
    #***********************************
    for x in range(width_roi):
        for y in range(height_roi):
            color = rgb_det[y % 2][x % 2]
            if color == 0: #R
                y_r = int(y/2)
                x_r = int(x/2)
                img_roi_r[y_r][x_r] = img_roi[y][x]
            elif color == 1: #Gr
                y_gr = int(y/2)
                x_gr = int((x-1)/2)
                img_roi_gr[y_gr][x_gr] = img_roi[y][x]
            elif color == 2: #Gb
                y_gb = int((y-1)/2)
                x_gb = int(x/2)
                img_roi_gb[y_gb][x_gb] = img_roi[y][x]
            else: #B
                y_gb = int((y-1)/2)
                x_gb = int((x-1)/2)
                img_roi_b[y_gb][x_gb] = img_roi[y][x]

    img_s[0][i] = np.mean(img_roi_r)
    img_s[1][i] = np.mean(img_roi_gr)
    img_s[2][i] = np.mean(img_roi_gb)
    img_s[3][i] = np.mean(img_roi_b)

    img_n[0][i] = np.std(img_roi_r)
    img_n[1][i] = np.std(img_roi_gr)
    img_n[2][i] = np.std(img_roi_gb)
    img_n[3][i] = np.std(img_roi_b)

    img_sn[0][i] = 20*(math.log10(img_s[0][i]/img_n[0][i]))
    img_sn[1][i] = 20*(math.log10(img_s[1][i]/img_n[1][i]))
    img_sn[2][i] = 20*(math.log10(img_s[2][i]/img_n[2][i]))
    img_sn[3][i] = 20*(math.log10(img_s[3][i]/img_n[3][i]))

fig1 = plt.figure()
ax1 = fig1.add_subplot(111)
ax1.plot(x_axis[0,:], img_sn[0,:],label = "R",color = "r")
ax1.plot(x_axis[0,:], img_sn[1,:],label = "Gr",color = "lightgreen")
ax1.plot(x_axis[0,:], img_sn[2,:],label = "Gb",color = "forestgreen")
ax1.plot(x_axis[0,:], img_sn[3,:],label = "B",color = "b")
ax1.set_xscale("log")
ax1.legend()

fig2 = plt.figure()
ax2 = fig2.add_subplot(111)
ax2.plot(x_axis[0,:], img_s[0,:],label = "R",color = "r")
ax2.plot(x_axis[0,:], img_s[1,:],label = "Gr",color = "lightgreen")
ax2.plot(x_axis[0,:], img_s[2,:],label = "Gb",color = "forestgreen")
ax2.plot(x_axis[0,:], img_s[3,:],label = "B",color = "b")
ax2.plot(x_axis[0,:], img_n[0,:],linestyle="--",label = "R_Noise",color = "r")
ax2.plot(x_axis[0,:], img_n[1,:],linestyle="--",label = "Gr_Noise",color = "lightgreen")
ax2.plot(x_axis[0,:], img_n[2,:],linestyle="--",label = "Gb_Noise",color = "forestgreen")
ax2.plot(x_axis[0,:], img_n[3,:],linestyle="--",label = "B_Noise",color = "b")
ax2.set_xscale("log")
ax2.set_yscale("log")
ax2.legend()

plt.show()