無料で人工衛星画像のタイムラプスの作り方．

１．はじめに

これまで，欧州の人工衛星であるSentinelの観測画像の入手する方法や，それをAPIを使って自動的に行う方法などを紹介してきました．

無料で最新の衛星画像を入手する方法．

人工衛星（Sentinel-2)の観測画像をAPIを使って自動取得してみた．

今回は，上記のコードをベースに，任意の場所の月スケールのタイムラプスのGifアニメーションを作ってみます．
例えば，２０１９年１月から１２月までの，Sentinel-２の観測画像より作成した東京付近の衛星画像のタイムラプスはこのようになります．

　（Credit：Europe Space Agency）　注：WEBに掲載するために画像サイズは落としています．

これを見ると，新国立競技場の屋根が徐々にできてきているのがわかります．

それでは，人工衛星Sentinel-２の観測画像によるタイムプラスの作成方法について紹介します．

コードはgitに置きましたので，個人の利用の範囲でよろしければご使用ください．
RS_Sentinel2_timelapse

２．　人工衛星画像のタイムプラスの作り方．

2.1 人工衛星画像の入手方法．

対象とする人工衛星は，欧州が開発・運用するSentinel-2になります．　
Sentinel-2の観測画像は，アカウント登録すればどなたでも無料で入手することができます．詳しくは，以下の記事を参考にしてください．

無料で最新の衛星画像を入手する方法．

自分の関心がある地域の画像を検索し，そのなかで適切な観測画像を選びダウロードし，そして加工する，を何度も行うことは手間になりますので，
Sentinelシリーズの観測画像のデータハブサイト（Copernicus hub)が提供するAPIを用いて，関心域の衛星画像を自動で入手します．

その方法，およびコードについては以下の記事を参考にしてください．

人工衛星（Sentinel-2)の観測画像をAPIを使って自動取得してみた．

上記のサイトのコードを加工することで，関心域のタイムプラス（一定時間毎の画像のアニメーション）を作ってみます．

２．２　人工衛星画像の選択．

人工衛星の観測画像の３割が雲が多く含まれているため，関心地域が雲に隠れて見えないことがあります．特に，みなさんも感じられていると思いますが，梅雨以降の秋口までは，台風などによって天気が曇天の日々が続くと，定期的に観測するSentinel-2の観測画像には雲が多く含まれます．
そのため，　人工衛星（Sentinel-2)の観測画像をAPIを使って自動取得してみた．で紹介したときと同じく，各月の観測画像の中から，一番雲のない画像を選択します．

products_gdf = api.to_geodataframe(products)
products_gdf_sorted = products_gdf.sort_values(['cloudcoverpercentage'], ascending=[True])
products_gdf_sorted

Sentinel-2の観測画像は１００km角で区切られており，ある程度オーバラップするようになっているため，画像の選択に注意が必要です．

前回の記事でも紹介しましたが，東京付近ですと千葉付近で提供画像がオーバラップしているのがわかります．
関心域をオーバラップするように選択すると，両方の画像がリストに入り，検索して雲がすくない画像でソートすると，関心域が少ない画像が選ばれることがあります．
その例が，前回の記事での紹介になります．（人工衛星（Sentinel-2)の観測画像をAPIを使って自動取得してみた．）

関心域を多く含む画像番号（タイル番号）を確認し，画像を取得することにします．
タイル番号は，以下のように設定されていますが，毎回確認することは手間になるため，取得した画像からタイル番号を取得し，それを設定することで対処します．

(Credit:European Space Agency)
Sentinel-2 product type

コードには説明を記載していますが，関心域で一番古い月の画像を確認用として取得し，そのjpeg画像を確認します．
その画像が思ったものとことなるときは，次に雲が少ない画像を選択し再度画像を取得し，加工して確認します．
その画像が適切であれば，タイル番号をコピーします．

上記でも説明がありますが，タイル番号は取得したファイル名を見るとわかります．

S2B_MSIL2A_20190103T013039_N0211_R074_T54SUE_20190103T033839

ファイル名が上記であれば，後半のT54SUEがタイル番号になります．
これをコピーし，placenumberにペーストします．

これでタイル番号が確定しましたので，同じタイル番号の各月の一番雲の少ない画像を選択し，jpeg画像を作成します．
この処理は，前回の記事と同じになります．

２．3　画像に観測日とクレジットを挿入．

人工衛星のGifアニメーションをつくるにあたり，その観測日を記載しないといつのものかわからなくなります．
また，使用する画像にはクレジットを入れる必要があります．その処理を行います．

フォントファイルがすでにあればそれを使用することでよいのですが，ここでは無料のフォントファイルを入手し利用させてもらいました．
いくつか無料でフォントを提供されていますが，今回は以下のサイトよりダウンロードします．（感謝です．）

M+ FONTS

#フォントファイルのダウンロード．
!wget https://osdn.net/dl/mplus-fonts/mplus-TESTFLIGHT-063a.tar.xz

#ダウンロードしたファイルの解凍．
!xz -dc mplus-TESTFLIGHT-*.tar.xz | tar xf -

#フォントファイルのアドレスを設定．
fontfile = "./mplus-TESTFLIGHT-063a/mplus-1c-bold.ttf"

ここで入手したフォントファイルをもちいて，観測日およびクレジットを以下のコードで挿入します．

#作成した衛星観測画像への撮像日とクレジットを記載
img = Image.open('./Image_jpeg_'+str(object_name) +'/' + str(Begin_date) + 'Masked_' +str(object_name) +'.jpg')

x = img.size[0]/100 #日付の記載位置の設定
y = img.size[1]/100 #日付の記載位置の設定
fs = img.size[0]/50 #日付のフォントサイズの設定
fs1 = int(fs)
obj_draw = ImageDraw.Draw(img)
obj_font = ImageFont.truetype("./font/mplus-1c-bold.ttf", fs1)
obj_draw.text((x, y), str(date), fill=(255, 255, 255), font=obj_font)
obj_draw.text((img.size[0]/2, img.size[1]-y - img.size[1]/20 ), 'produced from ESA remote sensing data', fill=(255, 255, 255), font=obj_font)

#作成した画像の保存．
img.save('./Image_jpeg_'+str(object_name) +'/' + str(Begin_date) + 'Masked_' +str(object_name) +'.jpg')

取得したSentinel-2の画像ファイルは１GBほどあるため，すべての画像を保存しているとストレージを圧迫するため，以下のコードで処理して削除しています．
（ダウンロードした画像ファイルを今後も利用したい場合は，以下のコードを削除してください．）

#Copernicus hubからのダウンロードファイル，および解凍フォルダの削除
    shutil.rmtree( str(product_title) + '.SAFE')
    os.remove(str(product_title) +'.zip')

これでGifアニメーションをつくるための画像の準備は終わりました．
次にアニメーションを作ります．

２．4　Gifアニメーションの作成．

先程作成したjpeg画像のアニメーションを作ります．画像のファイル名を画像の観測月名としたため，一番古い観測日にソートしてGifアニメーションを作成します．

images =[]

#画像ファイルのソート
files = sorted(glob.glob('./Image_jpeg_'+str(object_name) +'/*.jpg'))
images = list(map(lambda file: Image.open(file), files))

#gifアニメーションの作成と保存．
images[0].save('./Image_jpeg_'+str(object_name) +'/' + str(object_name) + '.gif', save_all=True, append_images=images[1:], duration=1000, loop=0)

これでgifの拡張子のファイルが作成されました．
変更時間（今回は１秒(duration =1000))などを変更するなど好みにあわせてください．
作成されたファイルをダウンロードし，ローカルにて動作を確認してください．

3. おわりに．

今回は,Sentinel-2の観測画像のタイムプラスを作成してみました．　開発が進んでいる地域・国であれば，より大きな変化がみられるかもしれません．
特に，湾岸地域は開発がすすんでいますので，そこを注目してみると面白いかもしれません．
経済指標との相関性はどうなのか，など考えてみると面白いかもしれません．

また，このコードをベースに差分をとってタイムスケールでの変化を抽出したり，使用するSentinel-2のバンドを組みあわせることで，植生の変化をみるのも面白いかもしれません．

みなさんの興味の参考になれば幸いです．　コメントなどいただければ嬉しいので，ご遠慮無くお願いします．
（タイル画像の設定する部分がスマートでないため，なんとかしたいです．．．）

今回も，全体のコードをgitに置きましたので，個人での利用の範囲で使ってください．
RS_Sentinel2_timelapse

参考記事

無料で最新の衛星画像を入手する方法．
人工衛星（Sentinel-2)の観測画像をAPIを使って自動取得してみた．
M+ FONTS