# -*- coding: utf-8 -*-
"""サハラの砂塵の動画（D1: 白紙・霞）。GEOS-FP の 3 時間ごとの場を、時間方向に線形補間して流す。止まりは置かない（連続する事象の時間）。"""
import sys, os, glob, shutil, subprocess
import numpy as np
sys.path.insert(0, '/home/claude/dust'); sys.path.insert(0, '/home/claude/forest')
import dust as Dm, forest1 as F
from PIL import Image
from scipy import ndimage
SS = 2; FPS_STEP = 8          # 3 時間を 8 コマ（1 日 ≈ 2.1 秒）
sc, lon, lat, ok, land = Dm.scene(SS)
PAPER = np.array([255,255,255],np.float32); LG = np.array([238,238,235],np.float32)
img0 = np.zeros((sc.H, sc.W, 3), np.float32); img0[:] = PAPER; img0[land] = LG
coast = (land ^ ndimage.binary_erosion(land, iterations=SS)).astype(np.float32)
F.NAVY = Dm.INK; F.WHITE = PAPER
base = np.array(F.comp_rel(Image.fromarray(img0.astype(np.uint8)), coast, 0.22, 0.16)).astype(np.float32)
ts = sorted(os.path.basename(p)[:-4] for p in glob.glob('/home/claude/dust/data/2025*.npz'))
print('steps', len(ts), ts[0], ts[-1])
def field(t): return F.sample_bilinear(ndimage.gaussian_filter(np.nan_to_num(Dm.grid('DUEXTTAU', t)), 1.0), lon, lat, ok)
FD = '/home/claude/dust/video/frames'; shutil.rmtree(FD, ignore_errors=True); os.makedirs(FD)
def compose(fr):
    al = (np.clip(fr/1.2, 0, 1)**0.8)[..., None]
    return Image.fromarray(np.clip(base*(1-al) + Dm.INK*al, 0, 255).astype(np.uint8)).resize(Dm.SIZE, Image.LANCZOS)
n = 0
f_prev = field(ts[0])
for _ in range(30): compose(f_prev).save(f'{FD}/{n:04d}.png'); n += 1          # 先頭 1 秒
for i in range(1, len(ts)):
    f_next = field(ts[i])
    for k in range(FPS_STEP):
        u = k / FPS_STEP; compose(f_prev*(1-u) + f_next*u).save(f'{FD}/{n:04d}.png'); n += 1
    f_prev = f_next
    if i % 8 == 0: print('day', i//8, n, flush=True)
last = compose(f_prev)
for _ in range(60): last.save(f'{FD}/{n:04d}.png'); n += 1
subprocess.run(['ffmpeg','-y','-loglevel','error','-framerate','30','-i',f'{FD}/%04d.png','-c:v','libx264','-pix_fmt','yuv420p','-crf','20',
                '-vf','pad=ceil(iw/2)*2:ceil(ih/2)*2:0:0:white','/home/claude/dust/video/saharan_dust_2025.mp4'], check=True)
print('frames', n, n/30, 's')