Raspberry Pi3で自動ノート取り装置を作った

はじめに

数理情報工学実験第二という演習で、Raspberry Piをつかって何かを作ることになりました。そこでAMATERASUという自動ノート取り装置を作ったので紹介します。

そもそもRaspberry Piって?

Raspberry Pi3 Model B ボード＆ケースセット (Element14版, Clear)-Physical Computing Lab

これです。安くて小型で色んなセンサーをつけて遊べるコンピュータです。今回はカメラモジュールを使いました。

自動ノート取り装置とは

自動ノート取りの目標は、講義を撮影した動画*1を処理することで、ノートの代わりとして使える画像を出力することです。具体的には次のgifのような画像を次々出力していくのを目標にしています。黒くなっている部分が何かについては後ほど説明しますが、そこには情報がないということを表しています。

動画を入力して処理することもできますが、オンラインで処理することを想定しています。

出力される画像は次のような条件を満たすようにします。

• 情報の漏れがないです。教員が黒板に書いたことは、全てどこかのタイミングで出力される画像に写っています。
• 無駄な情報を省きます。つまり、同じ板書が複数回出力されないようにします。
• 教員が写りこまないようにします。

基本戦略

AMATERASUがどのように出力する画像を作っていくかについて大まかに述べます。まず黒板を（だいたい）5秒おきに撮影します。その画像のうち、情報が新しくて、教員が映りこんでいない部分だけを残します。このような部分画像を次々足し合わせていくことで黒板全体が写った画像を得ます。また、黒板全体の画像が完成していて今出力すべきか否かについての判定も行います。

名前の由来

誓約で身の潔白を証明した建速須佐之男命は、高天原に居座った。そして、田の畔を壊して溝を埋めたり、御殿に糞を撒き散らしたりの乱暴を働いた。他の神は天照大神に苦情をいうが、天照大神は「考えがあってのことなのだ」とスサノヲをかばった。
しかし、天照大神が機屋で神に奉げる衣を織っていたとき、建速須佐之男命が機屋の屋根に穴を開けて、皮を剥いだ馬を落とし入れたため、驚いた1人の天の服織女は梭（ひ）が陰部に刺さって死んでしまった。ここで天照大神は見畏みて、天岩戸に引き篭った。高天原も葦原中国も闇となり、さまざまな禍（まが）が発生した。
天岩戸 - Wikipedia

「天岩戸に引きこもったアマテラスさんでもノートが手に入って卒業できるよ！」という理由です。あまり大した理由でもないし、頭文字をとるとAMATERASUになるようなかっこいい正式名称があるわけでもないです...

Raspberry Piを使う利点

Raspberry Piを使うことには以下のような利点があります。

• 安いです。カメラモジュールを含めても6000円～8000円程度で済みます。
• 小型なので比較的設置スペースの問題が起こりにくいです。
• ネットワークがなくてもその場で処理が行えます。
• ネットワークがあれば、dropboxã‚„twitterにアップロードでき、速報性があります。
• 今回はやっていませんが、他のセンサーを使えるという拡張性があります。

既存手法との比較

他の既存手法と比べて、どのようなメリットがあるのかを説明します。

å‹•ç”»

既存手法の一つ目は講義を録画することです。これと比べると、AMATERASUは

• 容量が少なくて済みます。2時間ほどの講義でも画像十数枚で済むので経済的です。
• 後で見直すときに便利です。動画より画像の方が簡単に講義全体を見渡せ、見たい部分をすぐに探せます。
• 本質的な情報は全て出力されています。確かに動画は情報を完璧に保存していますが、冗長でありAMATERASUで十分です。

パシャニキ

パシャニキとは近年大学で問題になっている黒板をスマホで撮影する人のことです。パシャニキと比べると、AMATERASUは

• シャッター音がしないので、他の受講者の迷惑になりません。突然立ち上がって視界をふさぐこともありません。
• 固定した状態で撮影されるので、複数の画像を比べても位置がそろっていてきれいです。手振れもありません。
• 教員が写りこまないので見やすいです。
• 撮影は自動で行われるため、講義を聴くことに集中できます。

ハードの解説

使ったもの

上で紹介したraspberry pi3を使っています。カメラモジュールは前年度のものを引き継いだため、具体的に何かは分からないです。ごめんなさい。

撮影方法

図のように黒板から数メートル離れた机の上にRaspberry piを設置します。その後カメラモジュールを割り箸とガムテープを使って固定します。撮影中机やカメラモジュールに衝撃が加わらないよう配慮されています。

一番初めには人も文字も写っていない背景画像が撮影されます。

ソフトの解説

環境

Python2を用いました。画像処理にはOpenCV 3.1を使いました。実行はiPython notebookを用いました。UbuntuとWindowsで動作確認しています。

概要

AMATERASUの概要について解説します。メモリ上にはoldとnowという2つの画像があります。oldは一回前に出力された画像で、nowは現在合成中の画像です。5秒おきにカメラモジュールが黒板の撮影し、capturedという画像が入ってきます。capturedは教員が映りこんでおらず、なおかつ新しい部分のみが切り出されます。その後、この部分画像がnowに足し合わされます。足し合わせるとき、すでにnowに含まれている領域だった場合はcapturedで上書きをします。また出力判定も行われ、必要ならHDDに保存したりツイートしたりします。

AMATERASUではnowや切り出されたcapturedのように部分画像を扱います。これらは1次元的にx座標（横方向）で区間分けして、有効かどうかを管理しています。同じx座標を持つピクセルは全て有効か全て無効かの2択です。2次元的に扱わないのは、計算を軽くするのと実装を簡単にするためです。

解決すべき問題

教員の除去

基本的には黒板を撮影すると教員が写りこんでしまいます。きれいな合成画像を得るためには、きちんと教員が写りこんでいる部分を識別して除去する必要があります。

出力判定・情報の新旧の判定

AMATERASUは現在出力すべきかどうかを判断する必要があります。また、出力が冗長にならないようにするには、古い情報を取り込まないよう気をつける必要があります。例を見ながらこのことを説明していきます。

黒板が3分割して使われている状況を考えます。このとき、教員ははじめに①、②、③という板書を左から順に行いました。その後、①を消して、そこに④という新しい板書をしました。②、③はまだ消されずに残っています。

このとき①を消して④を書き始めた段階で1枚目の板書は完成しています。その上、④をnowに取り込んでしまうと①が上書きされて消えてしまいます。そのためこのタイミングで出力を行わなくてはいけません。

また、④を書いたときに写真を撮ると、②、③という古い板書も一緒に写ってしまいます。もしこれがnowに足しあわされてしまうと、2枚目の出力画像にも②、③が含まれてしまい冗長になってしまいます。またこれは、後に紹介する出力判定のアルゴリズムを誤動作させる原因になります。そのため板書が古いものであることを判定し、もしそうならばnowに取り込まない仕組みが必要です。

このような処理は人の目から見れば当たり前なので簡単に思うかもしれませんが、Raspberry Piにやらせるためには機械的に扱える条件まで落とし込んでやる必要があります。

メインループ

AMATERASUは上に述べた問題を踏まえ画像を合成するために、5秒おきに新しい画像が撮影されるたび次のようなループを1回まわします。

画像が入ってくるとまず教員が検出され、その領域が消されます。次に板書が新しいか判定され、古いものは全て消されます。このようにして得られた部分画像を用い、出力判定が行われます。出力しないときはnowにこの部分画像が足し合わされます。

黒板領域の検出

ここから具体的な処理について述べていきます。まずは前処理についてです。撮影した画像は黒板以外のものも映っていますが、これらは必要ではありません。それどころか誤動作の原因にもなりえます。そのため黒板領域を検出してそれ以外をカットします。

これには一番初めに撮影した背景画像を用います。背景画像をHSV空間に変換し、緑色の領域が1になるよう2値化します。その後、1のピクセルが多いy座標の区間を黒板領域とします。

この黒板領域の検出はblackboardDetectionという関数により行われます。また、切り出す領域が確定した後は、撮影した画像はまずconvertToBlackboardSizeByValues関数によって黒板領域だけが切り出されます。

今回は撮影するときのカメラの固定の都合上画像が反転しているため、前処理で元の向きに戻しています。

教員検知

教員の検出について説明します。プログラム中ではfind_intervalという関数によって行われています。次の画像を例に説明します。個人情報保護のため赤色の四角で顔を隠しています。

まず画像を背景画像との絶対差分(ピクセルごとに差の絶対値を取る)を行い、その後2値化します。2値化する閾値は、画像を見て文字が書いてあるところがはっきりするよう定めました。*2

こうすると背景と異なる部分、つまり文字と人が白くなります。この状態から人と文字をうまく分離することができれば、人の検出ができます。人と文字との違いは、後者は細いことです。そこで縮小処理（erosion）を使うことで文字のみを消すことができます。縮小処理とは次の図のように、図形の外側を削り小さくします。

文字は細いため縮小処理によって消えてしまいますが、教員は太いため残ります。縮小処理をしたあと、拡大処理で元に戻すと次のような画像が得られます。

確かに人のみが残っています。黒板消しも残ってしまうため、教員検出の時は黒板の中央部のみを用います。その後、人がいないところのみ1を取るような配列を作ると画像の下のグラフが得られます。教員が写りこむのを避けるため、大きめに区間を取るようにしています。

画像差分による出力・更新判定

限定された減算

出力や更新の判定を行うとき、新しく書き込まれた場所や板書が消された場所といった情報が有用です。これらの情報を取るために、画像の差分を取ります。そのとき、絶対差分ではなく次のような「限定された減算」*3を用います。

引き算をした結果が負の数になるときは結果が0に丸められることだけが普通の減算と違います。画像に使う際はピクセルごとにこの演算を適応します。

絶対差分との比較

さて、前節で定義した減算は絶対差分と比べて何が利点なのでしょうか。

次のような例で考えてみます。黒板の左側に「A」と書かれている画像から、右側に「B」とかかれている画像を引きます。文字の部分が正になるよう2値化しているとします。このとき、絶対差分を取るとAとBの両方が残ります。それに対して、限定された減算を使うと、Aのみが残りBは残りません。Bの部分が引かれた結果は負になり0に丸められてしまうからです。

この非対称性により書き込みや削除の識別ができます。上のような引き算を考えます。左の画像の方が右のものより新しいとします。この状況ではAという新しい書き込みがあります。またCという板書が消されています。Bは変化していません。絶対差分を使うとAとC両方がでてきてしまいます。一方限定された減算を使うとAのみを取り出すことができます。つまり新しく書き込まれたところだけを取り出せています。左の画像のほうが古いときは時系列が逆になり、削除された情報Aのみを取り出すことができます。

実際の画像での例は、次節「出力判定への応用」で見ていきます。

プログラム中では、背景画像と絶対差分を取り2値化した後ピクセルごとに限定された減算を適応します。実際の画像は定点観測しているものの、微小な振動のせいで数ピクセル程度の位置のずれがあります。位置ずれのせいで演算の結果が正になると、書き込みや削除の誤検出につながります。これを避けるため、 の右側の画像には2値化後拡大処理を行っています。

出力判定への応用

出力すべきタイミングはいつなのか考えていきます。また黒板が三分割されている状況を例にします。まず①、②、③がかかれます。その後①を消して④を書きました。このとき１回目の出力をすべきです。その後②と③が消され⑤と⑥がかかれます。最後に④が消され⑦が書かれます。ここで2回めの出力をすべきです。

この例から、いつ出力すべきかを決める上で、板書が消されたという情報が有用だということが分かります。ただし、②と③が消された時には出力を行うべきではありません。まだ④しか書かれておらず、無駄な出力になってしまうからです。

以上をまとめると、出力すべきなのは「前回の出力後に書かれた部分が消されたとき」となります。前回の出力後に書かれた部分は、現在合成中のnowそのものです。削除の判定は、古い画像から新しい画像を限定された減算で引けばできました。つまり を計算して、正になるピクセルがあれば出力すべきとなります。この条件は特に左から書き始めることを仮定していません。なので、教員が左側の板書を残して新しい板書を始めた場合などでも、きちんと動作します。

実際の画像に用いると、この演算ではノイズが残ってしまいます。そこで、ガウスぼかしでノイズ除去を行っています。また古い情報の消去がうまくいかず間違えてnowに取り込んでいた場合、これが消されたとき誤出力をしてしまいます。そこでプログラム中では、「nowのうちoldの板書が消された部分」、つまり と共通部分を持つことも条件に加えることでロバストに動作するようにしています。

また教員は時に書き直しを行います。書き直しでも板書を消しているので は正になってしまいます。このときに出力すると無駄に画像の枚数が増えてしまいます。そこで、閾値をうまく設定することで1行程度消した場合には出力しないようにしています。

最後に実際の画像ではどうなるかを見てみましょう。1枚目がnowで2枚目がcapturedから教員を消去したものです。capturedの左側は教員が居るため、無効な区間として扱われ黒く表示されています。nowの黒い部分は、まだcapturedが足しあわされておらず、情報がなく無効な区間です。黒板の左側の板書がcapturedでは消されているのが分かります。3枚目の において、きちんとこの板書の消去が検出されています。最後にこれをx軸方向に射影して、閾値を越えたところを1とすることで4枚目の区間表現を得ます。

更新判定への応用

撮影した板書のうち、古い部分を消去して新しい部分のみをnowに足しあわせなくてはいけません。どう情報の新旧を判定しているかを見ていきます。

基本的には、capturedのうちoldと比べて新しく書き込まれた部分が新しい情報です。これは

で求めることができます。しかし、書き足しが行われると問題が生じます。次のような例を考えています。

oldでは①、②、③という板書がされていました。capturedでは①と②は消され、①の部分に④が書かれています。また、③の内容を④を書き始めた後付け足して、③'になっています。このような書き足しは時々起こります。このときwitternは③が書かれていた黒板の右側でも正になります。しかし、③は古い情報なのでnowに取り込むのは間違いです。oldの方に修正を加えるべき状況です。

そこで、「oldのうちcapturedでは消去された部分」を用います。これは

で求めることができます。消去された場所に書き込まれた情報は新しい情報だと考えられます。そこで、erasedとwrittenの共通部分が新しい情報だとします。このようにするとwrittenのうち、本当に新しい部分と書き足しを分けることができます。writtenが正となる区間ごとに、erasedと共通部分を持つか判定して、そうであるならばnowに加え、そうでないなら書き足しなのでoldに修正を加えます。

理想的には上のように書き足しをうまく識別して、oldを修正して再出力をすべきでした。実際には過度にoldの修正がおこりうまく働かなかったため、書き足しは無視して捨てるようにしています。この問題については後で詳しく述べます。

結果

nowの更新

AMATERASU～Automatic Blackboard Copy Device～
黒い部分は無効（まだ足されていない）部分。時々手とか紙とか残ってますけど最終的にうまく行ったし結果オーライで。

出力ファイル

時々黒い部分が残っていますが、そこはそもそも何も書き込まれていない部分なので問題ないです。今回は説明のため黒く残していますが、背景画像で黒い部分を埋めれば見栄えが良い画像が得られます。

問題点

試験回数の少なさ

今回は2回、別の人に対して同じ教室で実験を行いました。他の教室で正常に動くかは分かりません。また、他の人で試した場合、教員が今回想定していない動きをしたせいでうまく動かない可能性があります。*4

また、実験期間が足りなかったためオンライン処理の実験は行えていません。仕組み的におそらくオンラインでも動くとは思いますが...

環境依存性

いくつかの閾値を勝手に設定しているため、照明の具合によっては動かないかもしれません。またホワイトボードや移動式黒板には対応していません。

oldの更新問題

本当は書き足しがあったときは、oldを修正して再出力するようにしようとしていました。しかし、実際に動かしたところ、本来nowに入るべきものがoldの修正につかわれてしまい、oldがめちゃくちゃになるという問題が発生しました。そのため、この機能は削除されています。

原因はwrittenがきれいに求まるとは限らないことです。人間は画像を見てどこが固まりかを認識できますが、単純に射影して閾値を越えるかで判断しているAMATERASUではうまくいかない場合があります。平滑化をかけるようにはしていますが、本来ひとまとめであるべき部分がばらばらになったり、逆に離れているべき部分がくっついたりしてしまいます。前者が起こると、nowに足すべき画像がoldの修正に使われてしまいます。射影と閾値以外のうまい1次元化を使うべきでしょう。また、oldを修正するかどうかより厳しい基準で判断することでうまくいくかもしれません。

その他

1次元的に処理しているため、黒板を斜めに使われると誤動作してしまいます。これに関してはそういう書き方する教員が悪いと思いますが...

今回はカメラを固定していましたが、それでもブレが時折発生していました。せっかくいろいろなセンサーを使えるので加速度センサーでブレを検出すると面白いかなと思ってます。あとはずれたあとの位置あわせや、背景画像なしでの処理ができると、より汎用的になると思います。

まとめ

うまく動くものができた。すごい。

謝辞

@co_co_cocoa氏には、ゼミの撮影と実験への使用に快くご協力いただきました。ありがとうございます。
ココア (@co_co_cocoa) | Twitter
この実験は5人で行ったものであり、自分1人の成果ではありません。班員の皆さんにはアルゴリズムを考えるときの議論や実装時に大変お世話になりました。
数理情報第7研究室には、Raspberry Pi等の機材の購入、実験の指導、発表練習など様々な面で支援していただきました。ここに感謝の意を表します。

ソースコード

# -*- coding: utf-8 -*-
#get_ipython().magic(u'matplotlib inline')
#%matplotlib inline
import matplotlib.pyplot as plt
import numpy as np
import cv2
import time
import os
#from picamera.array import PiRGBArray
#from picamera import PiCamera

LECTURE_TIME=105*60
MODE=1 #0 camera 1 file
PICTURES_DIRECTORY="./Serial_Experiments_Cocoa/"
OUTPUT_DIRECTORY="./BlackboardRecord/"
INTERVAL_TIME=0
RECAPTURE_TIME=2
BLACKBOARD_WIDTH=1920

def convertToBlackboardSize(image):
    image2 = cv2.flip(image,-1)
    return image2[350:900,:,:]

def blackboardDetection(image):
    piet = image
    piet_hsv = cv2.cvtColor(piet, cv2.COLOR_BGR2HSV)
     
    # threshold for hue channel in blue range
    blue_min = np.array([70, 30, 70], np.uint8)
    blue_max = np.array([150, 255, 255], np.uint8)
    threshold_blue_img = cv2.inRange(piet_hsv, blue_min, blue_max)
     
    #threshold_blue_img = cv2.cvtColor(threshold_blue_img, cv2.COLOR_GRAY2RGB)
     
    #plt.imshow(threshold_blue_img)
    #plt.show()

    col_sum = np.sum(threshold_blue_img, axis=1)
    #plt.plot(col_sum)
    #plt.show()

    left = 540
    right = 540

    for i in range(540):
        k = 540 - i - 1
        if col_sum[k] > 1000:
            left = k
        else:
            break

    for i in range(540):
        k = 540 + i
        if col_sum[k] > 1000:
            right = k
        else:
            break

    return left, right

def convertToBlackboardSizeByValues(lower, upper, image):
    image2 = cv2.flip(image[lower:upper,:],-1)
    return image2
  
BACKGROUNDIMAGE_temp=cv2.imread(PICTURES_DIRECTORY+"blackboard0.jpg")
lower_cut, upper_cut = blackboardDetection(BACKGROUNDIMAGE_temp)
BACKGROUNDIMAGE = convertToBlackboardSizeByValues(lower_cut, upper_cut, BACKGROUNDIMAGE_temp)

# def find_interval(aaa):
#     kernel = np.zeros((4,4),np.uint8)
#     kernel[2,:]=1
#     kernel[3,:]=1
    
#     ddd7 = np.zeros(50*1920).reshape([50,1920])
#     ddd8 = cv2.threshold(aaa[100:495,:,0],140,188,cv2.THRESH_BINARY)[1]
#     ddd8 = cv2.erode(ddd8,kernel,iterations=3)
#     ddd8 = cv2.dilate(ddd8,kernel,iterations = 3)
    
#     ddd5 = np.concatenate([ddd7,ddd8])
#     hsv = cv2.cvtColor(aaa, cv2.COLOR_BGR2HSV)
    
#     ddd2 = cv2.threshold(hsv[:,:,1],70,np.max(hsv[:,:,0]),cv2.THRESH_BINARY)[1]
#     ddd1 = cv2.threshold(hsv[:,:,0],20,188,cv2.THRESH_BINARY)[1]
#     ddd1 = cv2.threshold(ddd1,20,188,cv2.THRESH_BINARY_INV)[1]
     
    
#     ddd3 = cv2.erode(ddd2[50:495,:],kernel,iterations=4)
#     ddd3 = cv2.dilate(ddd3,kernel,iterations = 4)
    
#     ddd4 = cv2.erode(ddd1[50:495,:],kernel,iterations=4)
#     ddd4 = cv2.dilate(ddd4,kernel,iterations = 4)
    
#     eee = np.mean( np.uint(ddd3) + np.uint(ddd4)+np.uint(ddd5) ,0)
#     left = 0
#     right = 0
#     for i,j in enumerate(eee):
#         if j>8:
#             left = i-200
#             if(left<=0):
#                 left = 0
#             if(i>=5):
#                 break
#     for i,j in enumerate(eee[::-1]):
#         if j>8:
#             right = 1920-i+200
#             if(right>1920):
#                 right = 1920
#             if(i>=5):
#                 break
#     ####print left,right
#     interval = np.ones(1920)
#     interval[left:right]=0
    
#     ddd = np.zeros(1920*550).reshape([550,1920])
#     ddd[:,:left]=0
#     ddd[:,right:1920]=0
#     ddd[:,left:right]=1
#     ddd = np.uint8(ddd)
#     ###true_image = cv2.bitwise_and(aa3,aa3,mask=(1-ddd))
#     #cv2.imwrite("./hoge/hey.jpg",aaa)
#     #plt.imshow(aaa)
#     return interval

def find_interval(image):
    diff = cv2.absdiff(image, BACKGROUNDIMAGE)
    grayed = cv2.cvtColor(diff, cv2.COLOR_BGR2GRAY)
    under_thresh = 30
    maxValue = 255
    th, drop_back = cv2.threshold(grayed, under_thresh, maxValue, cv2.THRESH_BINARY)

    eroded = cv2.erode(drop_back[50:450,:], neiborhood8, iterations = 3)
    dilated = cv2.dilate(eroded, neiborhood8, iterations = 15)

    row_sum = np.sum(dilated, axis=0)

    binarized = np.fmax(np.zeros(row_sum.size), thresold_1dim(row_sum,1000))
    res = np.ones(1920) - interval_expansion_strong(binarized, 300)
    return np.uint8(res)



def interval_expansion_strong(interval, iteration):
    width = iteration
    kernel = np.ones(iteration)
    expanded = np.convolve(interval, kernel, mode='same')
    res = np.fmax(np.zeros(interval.size), thresold_1dim(expanded, 1))
    return np.uint8(res)

SMOOTH_WIDTH=50
IMAGE_DIFFERENCE_TO_ARRAY_DILATION_ITERATION = 4
neiborhood8 = np.ones((3,3),dtype=np.uint8)

def filtered_erase(img1, img2):
    img = erase_detection(img1,img2)
    blured = cv2.GaussianBlur(img,(5,5),0)
    res = binarize(blured, 50)
    return res

def thresold_1dim(array, threshold):
    delta = array - threshold*np.ones(array.size)
    return np.sign(delta)

def binarize(img, threshold):
    under_thresh = threshold
    maxValue = 255
    th, drop_back = cv2.threshold(img, under_thresh, maxValue, cv2.THRESH_BINARY)
    return drop_back

def image_difference(img1, img2): #Compute max(img1 - img2, 0) and return 1 dimentional array. background is blackboard0.jpg
    diff1 = cv2.absdiff(img1, BACKGROUNDIMAGE)
    diff2 = cv2.absdiff(img2, BACKGROUNDIMAGE)
    erase = filtered_erase(diff1, diff2)
    return erase

def image_difference_to_array(img1, img2, interval,threshold): #Compute max(img1 - img2, 0) and return 1 dimentional array. background is blackboard0.jpg
    erase=image_difference(img1,img2)
    row_sum = np.sum(erase[50:450,:], axis=0)
    filtered = row_sum*interval
    kernel = np.ones(SMOOTH_WIDTH)/SMOOTH_WIDTH
    smoothed = np.convolve(filtered, kernel, mode='same')
    res = np.fmax(np.zeros(smoothed.size), thresold_1dim(smoothed,threshold))
    return np.uint8(res)

def strict_image_difference_to_array(img1, img2,interval, threshold):
    diff = image_difference(img1, img2)
    row_sum = np.sum(diff[50:450,:] ,axis=0)
    filtered = row_sum*interval
    binarized = np.fmax(np.zeros(filtered.size), thresold_1dim(filtered,threshold))
    smoothed = interval_expansion(interval_expansion(binarized, 40) , 40)
    return np.uint8(smoothed)

def image_difference3(img1, img2, img3): #Compute max(img1 - img2 - img3, 0) and return 1 dimentional array. background is blackboard0.jpg
    diff1 = cv2.absdiff(img1, BACKGROUNDIMAGE)
    diff2 = cv2.absdiff(img2, BACKGROUNDIMAGE)
    diff3 = cv2.absdiff(img3, BACKGROUNDIMAGE)

    grayed = cv2.cvtColor(diff2, cv2.COLOR_BGR2GRAY)
    under_thresh = 30
    maxValue = 255
    th, drop_back = cv2.threshold(grayed, under_thresh, maxValue, cv2.THRESH_BINARY)
    d_diff2 = cv2.dilate(drop_back,neiborhood8,iterations=IMAGE_DIFFERENCE_TO_ARRAY_DILATION_ITERATION)

    grayed = cv2.cvtColor(diff3, cv2.COLOR_BGR2GRAY)
    under_thresh = 30
    maxValue = 255
    th, drop_back = cv2.threshold(grayed, under_thresh, maxValue, cv2.THRESH_BINARY)
    d_diff3 = cv2.dilate(drop_back,neiborhood8,iterations=IMAGE_DIFFERENCE_TO_ARRAY_DILATION_ITERATION)

    grayed = cv2.cvtColor(diff1, cv2.COLOR_BGR2GRAY)
    under_thresh = 30
    maxValue = 255
    th, drop_back = cv2.threshold(grayed, under_thresh, maxValue, cv2.THRESH_BINARY)
    b_diff = drop_back

    delta = b_diff - d_diff2 - d_diff3
    res = (delta + abs(delta))/2

    blured = cv2.GaussianBlur(res,(5,5),0)
    res2 = binarize(blured, 50)
    return res2

def image_difference_to_array3(img1, img2, img3, interval,threshold): #Compute max(img1 - img2 - img3, 0) and return 1 dimentional array. background is blackboard0.jpg
    res2=image_difference3(img1,img2,img3)

    row_sum = np.sum(res2[50:450,:], axis=0)
    filtered = row_sum*interval
    kernel = np.ones(SMOOTH_WIDTH)/SMOOTH_WIDTH
    smoothed = np.convolve(filtered, kernel, mode='same')
    return np.uint8(np.fmax(np.zeros(smoothed.size), thresold_1dim(smoothed, threshold)))

def erase_detection(img1, img2): #TESTED ONLY FOR ABSOLUTE DIFFERENCE WITH BACKGROUND,detect lost of information of img1 when img2 is taken
    grayed = cv2.cvtColor(img2, cv2.COLOR_BGR2GRAY)
    under_thresh = 30
    maxValue = 255
    th, drop_back = cv2.threshold(grayed, under_thresh, maxValue, cv2.THRESH_BINARY)
    d_diff2 = cv2.dilate(drop_back,neiborhood8,iterations=IMAGE_DIFFERENCE_TO_ARRAY_DILATION_ITERATION)

    grayed = cv2.cvtColor(img1, cv2.COLOR_BGR2GRAY)
    under_thresh = 30
    maxValue = 255
    th, drop_back = cv2.threshold(grayed, under_thresh, maxValue, cv2.THRESH_BINARY)
    b_diff = drop_back

    delta = b_diff - d_diff2
    res = (delta + abs(delta))/2
    return res

def interval2intervals(interval):
    now1 = False
    a = -1
    intervals_pair = []
    for i in range(interval.size):
        if ((not now1) and interval[i] == 1):
            a = i
            now1 = True       
        elif (now1 and interval[i] == 0):
            intervals_pair.append((a,i))
            now1 = False      
        elif (now1 and i == interval.size - 1):
            intervals_pair.append((a,i+1))

    intervals = []
    for i in range(len(intervals_pair)):
        a = intervals_pair[i][0]
        b = intervals_pair[i][1]
        res = np.zeros(interval.size,np.uint8)
        for j in range(b-a):
            res[a+j] = 1
        intervals.append(np.array(res))

    return intervals
    
def interval_expansion(interval, iteration):
    kernel = np.concatenate((np.zeros(iteration),np.ones(iteration)))
    expanded1 = np.convolve(interval, kernel, mode='same')
    res1 = np.fmax(np.zeros(interval.size), thresold_1dim(expanded1, 1))
    kernel = np.concatenate((np.ones(iteration),np.zeros(iteration)))
    expanded2 = np.convolve(interval, kernel, mode='same')    
    res2 = np.fmax(np.zeros(interval.size), thresold_1dim(expanded2, 1))
    return np.uint8(np.fmax(interval,res1*res2))
    
def prodImgInterval(img,interval):
    interval2 =interval[:,np.newaxis]
    interval2 =interval2[np.newaxis,:]
    return img*np.uint8(interval2)


#mainLoop

#initialize
if MODE==0: 
    #camera=PiCamera()
    #rawCapture=PiRGBArray(camera)
    #time.sleep(0.1)
    starttime=time.time()
    nowtime=starttime
else:
    fileId=10
nowimage=np.array(BACKGROUNDIMAGE, copy=True)
nowInterval=np.zeros(BLACKBOARD_WIDTH,np.uint8)
oldimage=np.array(BACKGROUNDIMAGE, copy=True)
oldInterval=np.ones(BLACKBOARD_WIDTH,np.uint8)
firstOutputFlag=True
outputCount=-1
recaptureFlag=False
starttime=time.time()
nowtime=starttime
fileId=10
weakRenzoku=0

for f in os.listdir(OUTPUT_DIRECTORY[0:len(OUTPUT_DIRECTORY)-1]):
    os.remove(OUTPUT_DIRECTORY+f)
for f in os.listdir('./Debug'):
    os.remove('./Debug/'+f)

while nowtime<starttime+LECTURE_TIME:
    if MODE==0:
        #camera.capture(rawCapture,format="bgr")
        #image=rawCapture.array
        print("no camera")
    else:
        if fileId>8799:
            break
        image=cv2.imread(PICTURES_DIRECTORY+"blackboard%d.jpg"%fileId)
        fileId+=10
    
    image=convertToBlackboardSizeByValues(lower_cut, upper_cut,image)
    nohumanInterval=find_interval(image)
    #if recaptureFlag:
    #    if MODE==0:
    #        time.sleep(RECAPTURE_TIME)
    #    continue
    #oldgraph=prodImgInterval(oldimage,oldInterval*nohumanInterval)
    #imagegraph=prodImgInterval(image,nohumanInterval)

    #characteristic functions
    old_and_now = oldInterval*nowInterval
    only_old = oldInterval*(1-nowInterval)
    only_now = nowInterval*(1-oldInterval)
    neither_old_nor_now = (1-nowInterval)*(1-oldInterval)

    erasedInterval=image_difference_to_array(oldimage, image, oldInterval*nohumanInterval,200)*oldInterval*nohumanInterval
    old_minus_now_Interval = image_difference_to_array(oldimage, nowimage, nowInterval*oldInterval,200)*oldInterval*nowInterval
    #writtenInterval= (old_and_now*image_difference_to_array3(image, oldimage,nowimage,nohumanInterval*old_and_now) + only_now*image_difference_to_array(image,nowimage,nohumanInterval*only_old) + only_old*image_difference_to_array(image, oldimage,nohumanInterval*only_old) + neither_old_nor_now*image_difference_to_array(image, BACKGROUNDIMAGE, nohumanInterval*neither_old_nor_now))*nohumanInterval                        
    owrittenInterval = nohumanInterval*(oldInterval*image_difference_to_array(image,oldimage,nohumanInterval*oldInterval,200) \
                        + (1-oldInterval)*image_difference_to_array(image, BACKGROUNDIMAGE, nohumanInterval,200))
    writtenInterval = interval_expansion(owrittenInterval, 30)
    writtenInterval = interval_expansion(writtenInterval, 30)    
    writtenIntervals=interval2intervals(writtenInterval)
        
    if outputCount==-1:
        notCompatibleInterval=strict_image_difference_to_array(nowimage,image,nohumanInterval*nowInterval,3000)*nowInterval*nohumanInterval
        weakNotCompatibleInterval=image_difference_to_array(nowimage,image,nohumanInterval*nowInterval,200)*nowInterval*nohumanInterval
    else:
        notCompatibleInterval=strict_image_difference_to_array(nowimage,image,nohumanInterval*nowInterval,3000)*nowInterval*nohumanInterval*old_minus_now_Interval
        weakNotCompatibleInterval=image_difference_to_array(nowimage,image,nohumanInterval*nowInterval,200)*nowInterval*nohumanInterval*old_minus_now_Interval    
    s=np.sum(notCompatibleInterval)
    t=np.sum(weakNotCompatibleInterval)
    if s>100:
        #debug_output_all_status() 
        weakRenzoku=0
        if outputCount==-1:
            outputCount=0
        else:
            oldimage=prodImgInterval(oldimage,oldInterval)
            cv2.imwrite(OUTPUT_DIRECTORY+'output%d_%d.png'%(outputCount, fileId),oldimage)
            outputCount+=1
#            print('output%d fileId=%d'%(outputCount,fileId))            
#            plt.imshow(prodImgInterval(oldimage,oldInterval))
#            plt.show()
        oldimage=np.uint8(np.array(nowimage,copy=True))
        nowimage=np.uint8(np.array(BACKGROUNDIMAGE,copy=True))
        oldInterval = np.array(nowInterval,copy=True)
        nowInterval = np.zeros(BLACKBOARD_WIDTH)
    elif (weakRenzoku==0 and t>100) or (0<weakRenzoku and weakRenzoku<3):
        weakRenzoku+=1
    else:
        weakRenzoku=0
        nowWrittenInterval=np.zeros(BLACKBOARD_WIDTH)
        for writtenInt in writtenIntervals:    
            if np.sum(writtenInt*erasedInterval)>20 or outputCount==-1:
                #debug
                #DEBUG_OUTPUT_WRITTENINTERVAL_N_TIMES_UPDATE_COUNTER += 1
                #if DEBUG_OUTPUT_WRITTENINTERVAL_N_TIMES_UPDATE_LOWER <= DEBUG_OUTPUT_WRITTENINTERVAL_N_TIMES_UPDATE_COUNTER and DEBUG_OUTPUT_WRITTENINTERVAL_N_TIMES_UPDATE_COUNTER <= DEBUG_OUTPUT_WRITTENINTERVAL_N_TIMES_UPDATE_UPPER:
                    #debug_output_all_status()
                    #print('WRITTENINTERVAL was outputed after'+str(DEBUG_OUTPUT_WRITTENINTERVAL_N_TIMES_UPDATE_COUNTER)+'times update')
                #debug end
                nowimage=prodImgInterval(image,writtenInt)+prodImgInterval(nowimage,nowInterval*(1-writtenInt))
                nowInterval=1-(1-nowInterval)*(1-writtenInt)
                nowWrittenInterval+=writtenInt
            #else:
            #    oldimage=prodImgInterval(image,writtenInt)+prodImgInterval(oldimage,oldInterval*(1-writtenInt))
            #    oldInterval=1-(1-oldInterval)*(1-writtenInt)
            
    # print("fileId=%d s=%d oC=%d wR=%d"%(fileId,s,outputCount,weakRenzoku))
    # if fileId<1500:
    #     continue
    # Row=3
    # Column=3
    # fig,ax=plt.subplots(Row,Column,figsize=(10,4))
    # fig.tight_layout();
    # for r in range(Row):
    #     for c in range(Column):
    #         ax[r,c].axis([0,BLACKBOARD_WIDTH,-0.2,1.2])
            
    # ax[0,0].plot(notCompatibleInterval)
    # ax[0,0].set_title('notCompatibleInterval')
    
    # ax[1,0].plot(writtenInterval)
    # ax[1,0].set_title('writtenInterval')
    
    # ax[2,0].plot(nowWrittenInterval)
    # ax[2,0].set_title('nowWrittenInterval')
    
    # ax[0,1].plot(erasedInterval)
    # ax[0,1].set_title('erasedInterval') 
    
    # ax[1,1].plot(nohumanInterval)
    # ax[1,1].set_title('nohumanInterval')    
    # #ax[1,1].plot(nowInterval)
    # #ax[1,1].set_title('nowInterval') 

    # ax[2,1].plot(old_minus_now_Interval)
    # ax[2,1].set_title('old_minus_now_Interval') 
   
    # ax[0,2].plot(nowInterval)
    # ax[0,2].set_title('nowInterval') 
    
    # ax[1,2].plot(weakNotCompatibleInterval)
    # ax[1,2].set_title('weakNotCompatibleInterval')    

    # ax[2,2].plot(oldInterval)
    # ax[2,2].set_title('oldInterval')    
   
    #plt.show()
    plt.savefig('./Debug/intervals%d.png'%fileId)
    temptime=nowtime
    nowtime=time.time()
    sleeptime=INTERVAL_TIME-(nowtime-temptime)
    if sleeptime>0:# and MODE==0:
        time.sleep(sleeptime)

#last output
cv2.imwrite(OUTPUT_DIRECTORY+'output%d.png'%outputCount,prodImgInterval(oldimage, oldInterval))
outputCount += 1
cv2.imwrite(OUTPUT_DIRECTORY+'output%d.png'%outputCount,prodImgInterval(nowimage, nowInterval))

#mainLoopEnd

#debug
    
DEBUG_OUTPUT_WRITTENINTERVAL_N_TIMES_UPDATE_LOWER = 0
DEBUG_OUTPUT_WRITTENINTERVAL_N_TIMES_UPDATE_UPPER = 0
DEBUG_OUTPUT_WRITTENINTERVAL_N_TIMES_UPDATE_COUNTER = 0

def debug_output_all_status():
    print('fileID,s')
    print fileId,s
    print('oldimage')
    plt.imshow(prodImgInterval(oldimage, oldInterval))
    cv2.imwrite('./Debug/'+ str(fileId) + 'oldimage.png',prodImgInterval(oldimage, oldInterval))
    plt.show()
    print('nowimage')
    plt.imshow(prodImgInterval(nowimage,nowInterval))
    plt.show()
    cv2.imwrite('./Debug/'+ str(fileId) + 'nowimage.png',prodImgInterval(nowimage, nowInterval))
    print('image')
    plt.imshow(image)
    plt.show()
    cv2.imwrite('./Debug/'+ str(fileId) + 'image.png',image)
    erase=image_difference(nowimage,image)
    res=image_difference_to_array(nowimage,image,nohumanInterval*nowInterval)
    #temptemp=np.array(erase,copy=True)
    print('now_erase')
    plt.imshow(erase)
    cv2.imwrite('./Debug/'+ str(fileId) + 'now_erase.png',erase)
    plt.show()
    print('now_erased_non_binarized')
    row_sum = np.sum(erase[50:500,:], axis=0)
    width = 50
    filtered = row_sum*nohumanInterval*nowInterval
    kernel = np.ones(width)/width
    smoothed = np.convolve(filtered, kernel, mode='same')
    plt.plot(smoothed)
    plt.show()
    print('now_erased(Interval)')
    plt.plot(res)
    plt.show()
    print('nohumanInterval')
    plt.plot(nohumanInterval)
    plt.show()
    print('writtenInterval_image')
    writtenInterval_image = nohumanInterval*(oldInterval*image_difference(image,oldimage) + (1-oldInterval)*image_difference(image, BACKGROUNDIMAGE))
    #writtenInterval_image= old_and_now*image_difference_to_array3_(image, oldimage,nowimage) + only_now*image_difference_to_array_(image,nowimage) + only_old*image_difference_to_array_(image, oldimage) + neither_old_nor_now*image_difference_to_array_(image, BACKGROUNDIMAGE)
    plt.imshow(writtenInterval_image)
    cv2.imwrite('./Debug/'+str(fileId)+'writtenInterval.png', writtenInterval_image)
    plt.show()
    print('writtenInterval')
    plt.plot(writtenInterval)
    plt.show()

    erase=image_difference(oldimage,image)
    res=image_difference_to_array(oldimage,image,nohumanInterval*oldInterval)
    #temptemp=np.array(erase,copy=True)
    print('erased(old - image)')
    plt.imshow(erase)
    cv2.imwrite('./Debug/'+ str(fileId) + 'erase(old-image).png',erase)
    plt.show()
    print('erased_non_binarized')
    row_sum = np.sum(erase[50:500,:], axis=0)
    width = 50
    filtered = row_sum*nohumanInterval*oldInterval
    kernel = np.ones(width)/width
    smoothed = np.convolve(filtered, kernel, mode='same')
    plt.plot(smoothed)
    plt.show()
    print('erasedInterval')
    plt.plot(res)
    plt.show()

    erase=image_difference(oldimage,nowimage)
    res=image_difference_to_array(oldimage,nowimage,nowInterval*oldInterval)
    #temptemp=np.array(erase,copy=True)
    print('old - now')
    plt.imshow(erase)
    cv2.imwrite('./Debug/'+ str(fileId) + 'old-now.png',erase)
    plt.show()
    print('old_minus_now_non_binarized')
    row_sum = np.sum(erase[50:500,:], axis=0)
    width = 50
    filtered = row_sum*nohumanInterval*oldInterval
    kernel = np.ones(width)/width
    smoothed = np.convolve(filtered, kernel, mode='same')
    plt.plot(smoothed)
    plt.show()
    print('old_minus_now_Interval')
    plt.plot(res)
    plt.show()
#debug end

*1:正確には連続撮影された画像

*2:応用を考えるなら閾値を設定するのも自動化するか、この値で環境が変わってもうまくいくか確認すべきです。

*3:Wikipedia先生はこう呼んでました。計算理論とかでよくでてくるやつですね。 原始再帰関数 - Wikipedia

*4:今回でも人が持っている紙が検出しにくいとか、突然しゃがまれるといった予想していなかった問題がありました。