1 | #!/usr/bin/env python3 |
---|
2 | |
---|
3 | import matplotlib |
---|
4 | import sys |
---|
5 | from matplotlib import colors, ticker, cm |
---|
6 | from matplotlib.backends.backend_pdf import PdfPages |
---|
7 | from netCDF4 import Dataset |
---|
8 | import matplotlib.pyplot as plt |
---|
9 | import numpy as np |
---|
10 | import math |
---|
11 | from matplotlib import animation |
---|
12 | import matplotlib.gridspec as gridspec |
---|
13 | from matplotlib.lines import Line2D |
---|
14 | from matplotlib.collections import PatchCollection |
---|
15 | from matplotlib.patches import Polygon |
---|
16 | matplotlib.rcParams['text.usetex'] = True |
---|
17 | |
---|
18 | |
---|
19 | input_file_name = './canyon_agt.001.nc' |
---|
20 | |
---|
21 | file_content = Dataset(input_file_name, mode='r') # load entire netCDF file |
---|
22 | print(file_content) |
---|
23 | time = file_content.variables['time'][:] # extract variable |
---|
24 | agx = file_content.variables['ag_x'][:,:] |
---|
25 | agy = file_content.variables['ag_y'][:,:] |
---|
26 | grp = file_content.variables['ag_group'][:,:] |
---|
27 | file_content.close() # close netCDF file |
---|
28 | |
---|
29 | input_file_name = './topo.txt' |
---|
30 | f = open(input_file_name, mode='r') |
---|
31 | data = [] |
---|
32 | building_id = [] |
---|
33 | vertex_id = [] |
---|
34 | polygon = [] |
---|
35 | gons = [] |
---|
36 | p_last = 1 |
---|
37 | for line in f: |
---|
38 | line = line.strip() |
---|
39 | columns = line.split() |
---|
40 | if (int(columns[0]) > p_last): |
---|
41 | gons.append(polygon) |
---|
42 | polygon = [] |
---|
43 | polygon.append([float(columns[2]),float(columns[3])]) |
---|
44 | p_last = int(columns[0]) |
---|
45 | gons.append(polygon) |
---|
46 | nop = int(columns[0]) |
---|
47 | patches = [] |
---|
48 | |
---|
49 | for i in range(nop): |
---|
50 | poly = Polygon(gons[i-1], True) |
---|
51 | patches.append(poly) |
---|
52 | p = PatchCollection(patches,facecolors='grey') |
---|
53 | |
---|
54 | Nt = np.shape(time)[0] |
---|
55 | print(Nt) |
---|
56 | Nt = Nt -1 |
---|
57 | Nags = np.shape(agx[1,:]) |
---|
58 | Nx = 40 |
---|
59 | Ny = 40 |
---|
60 | |
---|
61 | x = np.linspace(0.0, 40.0, Nx) |
---|
62 | y = np.linspace(0.0, 40.0, Ny) |
---|
63 | X, Y = np.meshgrid(x, y) |
---|
64 | |
---|
65 | fig,ax = plt.subplots() |
---|
66 | fig.set_size_inches(4., 8.) |
---|
67 | xrange = [10, 30] |
---|
68 | yrange = [0, 40] |
---|
69 | |
---|
70 | |
---|
71 | def animate(i): |
---|
72 | ax.clear() |
---|
73 | ax.add_collection(p) |
---|
74 | x = agx[i,:] |
---|
75 | y = agy[i,:] |
---|
76 | ax.set_xlim(*xrange) |
---|
77 | ax.set_ylim(*yrange) |
---|
78 | ax.scatter(x,y,s=15,c=grp[i,:]) |
---|
79 | plt.title(r'simulated time: %i s' % time[i]) |
---|
80 | k = i/Nt*100. |
---|
81 | print("Done: %6.2f percent\r" % k, file=sys.stdout, end=" ") |
---|
82 | return ax, |
---|
83 | |
---|
84 | anim = animation.FuncAnimation(fig, animate, frames=Nt, interval=60, blit=False,save_count=0) |
---|
85 | anim.save('ags.mp4', dpi=200, writer='ffmpeg') |
---|