flow.py

import numpy as np
import os
import matplotlib.pyplot as plt
import matplotlib.cm as cm
from matplotlib.colors import LightSource


def image_dem(filename,gsx,gsy,pix,surface_dem):
    dpi = 300.0 #72.0
    # Create surface dem map
    solar_angle = 20.0
    dem_map = np.copy(surface_dem) - np.roll(surface_dem, 1, 0)
    dem_map = (0.5 * np.pi) + np.arctan2(dem_map, pix)
    dem_map = dem_map - np.radians(solar_angle) * (0.5 * np.pi)
    np.place(dem_map, dem_map > (0.5 * np.pi), 0.5 * np.pi)
    dem_map = np.absolute(dem_map)
    dem_map = 254.0 * np.cos(dem_map)

    # Save image to file
    height = gsy/ dpi
    width = gsx / dpi
    fig = plt.figure(figsize=(width, height), dpi=dpi)
    fig.figimage(dem_map, cmap=cm.gray, origin='lower')
    plt.savefig(filename)
    print(f'Image saved to {filename}')

    return

os.system(f'gfortran -o flow -O3 -fopenmp module_globals.f90 module_realistic.f90 realistic_perlin_noise.f90 flow.f90')

zoom = 1.0
xoff = 0.0
yoff = 0.0
damp_scale = 1.0
damp0 = 1.0
damp = 0.8
warp0 = 0.0
warp = 0.00
dem_file = 'flow.dat'
pix = 1.000
gsx = 1920
gsy = 1080
nframes = 2400
pdistort_scale = 10.0
for i in range(nframes):
    dist_pers = 1.5 * (i / nframes)
    pfile = open('flow.in', 'w')
    print(f'{gsx} {gsy}', file=pfile)
    print(f'{xoff} {yoff}', file=pfile)
    print(f'{i}',file=pfile)
    pfile.close()
    os.system(f'./flow')
    surface_dem = np.fromfile(dem_file , dtype = np.float64)
    surface_dem = surface_dem.reshape((gsy,gsx) , order='C')
    imgfile = f'flowframes/flow{i:06d}.png'
    image_dem(imgfile,gsx,gsy,pix,surface_dem)

os.system(f'open {imgfile}')