Scikit-learn 青空を再現描画

#始めに

画像などで、良く電線などが写っていて、何とかしたいと思うことがあります。
上記画像で電線など余計なものを、機械学習などの応用で消すことはできないか？
と思い、この鉄塔の右側を雲一つない晴天に書き換えることをScikit-learnを使って
やってみたいと思います。
Scikit-learnを使って技術的なお遊びで、気楽にやっていますので、scoreを取ったりとか
正解率がどうのということは記事中でははぶきます。

#処理
##処理概要
処理としては下記のように行いました。

1. 画像読み込み
2. 学習用データを作成
3. Pipelineで回帰計算インスタンスをRGBそれぞれ定義し学習させる。
4. 回帰により推定値をR,G<Bそれぞれ取得
5. 推定値で青空を描き、マスク処理で画像合成

##回帰計算
空のR,G,B値を求めて、空を描き直します。
機械学習として、座標x, 座標yを入力し、R,G,Bを出力するように
したいですが求められる推定値は1つなので、R,G,B毎に1つづつ
用意します。

• B ← f(座標x, 座標y)
• G ← f(座標x, 座標y)
• R ← f(座標x, 座標y)

**f()はPipeline()**でまとめられた回帰計算処理です。

##サンプリング
サンプリングは下記位置23か所で行います。

Pythonコード

mes_pnts = np.array([
    (174,  6),( 11, 12),(362, 24),(233, 28),( 66, 29),
    (207, 32),(  9, 54),( 96, 62),(342, 73),(289, 92),
    (375,112),(387,152),(172,155),(230,163),(191,168),
    (176,188),(345,195),(  9,222),(383,230),(291,236),
    (208,256),( 47,259),(156,262),])

# 機械学習用データを作成
rw,rh = 3, 3
x_r,x_g,x_b = [],[],[]
for ps in mes_pnts:
    x = int(ps[0]-1)
    y = int(ps[1]-1)

    # サンプリング座標を中心にして
    # 3x3の9画素の平均を取得する
    m = cv2.mean(img[y:y+rh, x:x+rw])
    x_r.append(m[2])    # R平均
    x_g.append(m[1])    # G平均
    x_b.append(m[0])    # B平均

サンプリング値は指定座標を中心に3×3の領域の平均値を使用しています。
画像処理で1画素で判断するのは危険です。

##Pipelineの使用
Pipelineは下記機能を組み合わせて使用します。

• StandardScalerの標準化機能
• PolynomialFeaturesの多項式機能
• LinearRegressionの回帰分析機能

これを使うことで入力、学習の操作がいっぺんにできます。
上記１つ１つで１記事が書けるほどなので説明がしきれず、本記事の守備範囲を
超えてしまうと感じる為、詳細についてはインターネット、Qiita上に沢山あるので
詳しく知りたい方はググって見て下さい。

参考
  scikit-learn.org - sklearn.pipeline.Pipeline
  scikit-learn.org - sklearn.preprocessing.PolynomialFeatures
  scikit-learn.org - sklearn.linear_model.LinearRegression
  scikit-learn.org - sklearn.preprocessing.StandardScaler


Rチャンネルを処理するPythonコード

from sklearn.pipeline import Pipeline
from sklearn.linear_model import LinearRegression 
from sklearn.preprocessing import PolynomialFeatures
from sklearn.preprocessing import StandardScaler

order = 5
pf_r = Pipeline([
        ("scaler", StandardScaler()),
        ("poly_features", PolynomialFeatures(degree=order,include_bias=False)),
        ("linear_reg", LinearRegression(fit_intercept=True)) ])
pf_r.fit(mes_pnts, x_r) # R学習

# 推定
prod_r = np.clip(pf_r.predict(wp),0,255).astype(np.int32)

##マスク処理
マスク画像としては下記画像を使います。白い部分(255)に青空を描きます。

マスク処理の流れ

1. 反転マスク画像 ← NOT マスク画像
2. 合成用画像 ← 入力画像 AND 反転マスク画像
3. 結果画像 ← 生成青空画像 AND 合成用画像

#Pythonコード

from sklearn.pipeline import Pipeline
from sklearn.linear_model import LinearRegression 
from sklearn.preprocessing import PolynomialFeatures
from sklearn.preprocessing import StandardScaler
import numpy as np
import math
import cv2

def main():
    # 鉄塔画像の読み込み
    img = cv2.imread('sky.png', cv2.IMREAD_COLOR)
    # マスク画像の読み込み
    mask = cv2.imread('rchmask.png', cv2.IMREAD_GRAYSCALE)
    mask = cv2.bitwise_not(mask)
    height, width = img.shape[:2]
    mes_pnts = np.array([   # サンプリング座標
        (174,  6),( 11, 12),(362, 24),(233, 28),( 66, 29),
        (207, 32),(  9, 54),( 96, 62),(342, 73),(289, 92),
        (375,112),(387,152),(172,155),(230,163),(191,168),
        (176,188),(345,195),(  9,222),(383,230),(291,236),
        (208,256),( 47,259),(156,262),])

    # 機械学習用データを作成
    rw,rh = 3, 3
    x_r,x_g,x_b = [],[],[]
    for ps in mes_pnts:
        x = int(ps[0]-1)
        y = int(ps[1]-1)
        # サンプリング座標を中心にして
        # 3x3の9画素の平均を取得する
        m = cv2.mean(img[y:y+rh, x:x+rw])
        x_r.append(m[2])    # R平均
        x_g.append(m[1])    # G平均
        x_b.append(m[0])    # B平均

    order = 5
    pf_r = Pipeline([
            ("scaler", StandardScaler()),
            ("poly_features", PolynomialFeatures(degree=order,include_bias=False)),
            ("linear_reg", LinearRegression(fit_intercept=True)) ])
    pf_r.fit(mes_pnts, x_r) # R学習
    pf_g = Pipeline([
            ("scaler", StandardScaler()),
            ("poly_features", PolynomialFeatures(degree=order,include_bias=False)),
            ("linear_reg", LinearRegression(fit_intercept=True)) ])
    pf_g.fit(mes_pnts, x_g) # G学習
    pf_b = Pipeline([
            ("scaler", StandardScaler()),
            ("poly_features", PolynomialFeatures(degree=order,include_bias=False)),
            ("linear_reg", LinearRegression(fit_intercept=True)) ])
    pf_b.fit(mes_pnts, x_b) # B学習

    # マスクより空を再描画する
    ya, xa = np.where(mask == 0)
    wp = np.hstack((xa.reshape(-1,1), ya.reshape(-1,1)))

    # 推定
    prod_r = np.clip(pf_r.predict(wp),0,255).astype(np.int32)
    prod_g = np.clip(pf_g.predict(wp),0,255).astype(np.int32)
    prod_b = np.clip(pf_b.predict(wp),0,255).astype(np.int32)

    # 結果画像
    result = cv2.bitwise_and(img, cv2.merge((mask,mask,mask)))
    result[ya, xa, 0] = prod_b
    result[ya, xa, 1] = prod_g
    result[ya, xa, 2] = prod_r

    # 表示
    cv2.imshow("original", img)
    cv2.imshow("result", result)
    cv2.waitKey(0)

if __name__ == '__main__':
    main()

sky.pngファイルが入力画像です。

#結果
上記コードを実行することで、下記画像が得られます。

#おわりに
画像としては、いいかんじにできたかなとは思います。際の部分ではイマイチ感ありますが。
この処理では、実は入力のサンプリング位置の確定や、マスク画像を作る方が
難しかったりします。その辺の方法についても、今後いろいろ試したいなとは思います。
お読み下さりありがとうございました。

#参考
Emotion Explorer - 機械学習で空いっぱいの電線、鉄塔を消すことを試す
Emotion Explorer - Scikit-learn Pipelineを試す
Emotion Explorer - Scikit-learn 多項式回帰と交互作用項

#写真引用
鉄塔と青空 photoAC掲載 ボッチ写真家さんの写真