Skip to content

Commit 1d44880

Browse files
mcgrattamcgratta
authored
Merge pull request #14821 from mcgratta/master
Python conversions
2 parents 958d897 + 30a5f99 commit 1d44880

File tree

3 files changed

+472
-0
lines changed

3 files changed

+472
-0
lines changed

‎Utilities/Python/scripts/FHWA_Tunnel.py‎

Lines changed: 198 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,198 @@
1+
2+
import numpy as np
3+
import pandas as pd
4+
import matplotlib.pyplot as plt
5+
from scipy.interpolate import griddata
6+
import os
7+
import warnings
8+
# include FDS plot styles, etc.
9+
import fdsplotlib
10+
11+
warnings.filterwarnings('ignore')
12+
13+
# McGrattan
14+
# 9-27-2022
15+
# FHWA_Tunnel.py
16+
#
17+
# This script creates several different kinds of contour and scatter plots for the FHWA Tunnel simulations.
18+
19+
# clear all - not needed in Python
20+
# close all
21+
plt.close('all')
22+
23+
outdir = '../../../out/FHWA_Tunnel/'
24+
expdir = '../../../exp/FHWA_Tunnel/'
25+
pltdir = '../../Manuals/FDS_Validation_Guide/SCRIPT_FIGURES/FHWA_Tunnel/'
26+
27+
pos = [1.0, 1.5, 2.0, 2.5, 3.0, 3.5, 4.0, 4.5, 5.0, 5.5, 6.0, 6.5, 7.5, 9.5, 10.5, 11.5]
28+
test = ['IFAB-07', 'IFAB-08', 'IFAB-09', 'IFAB-10', 'IFAB-11', 'IFAB-13', 'IFAB-14', 'IFAB-15', 'IFAB-19', 'IFAB-22', 'IFAB-24']
29+
test2 = ['Test 7', 'Test 8', 'Test 9', 'Test 10', 'Test 11', 'Test 13', 'Test 14', 'Test 15', 'Test 19', 'Test 22', 'Test 24']
30+
single_level = [50]
31+
setpoint = [10000, 400, 399, 338, 240, 322, 390, 420, 360, 10000, 10000]
32+
33+
plot_style = fdsplotlib.get_plot_style("fds")
34+
35+
def addverstr(ax, git_filename, style, x, y, fontname, interpreter, fontsize):
36+
"""Add version string to plot - placeholder implementation"""
37+
# This function would read git information and add it to the plot
38+
# For now, we'll add a placeholder text
39+
if os.path.exists(git_filename):
40+
ax.text(x, y, 'Git Version Info', transform=ax.transAxes,
41+
fontname=fontname, fontsize=fontsize)
42+
43+
for k in range(11): # Experiments (0-based indexing in Python)
44+
45+
n_res = 1
46+
# if k==7: # k==8 in MATLAB (1-based) corresponds to k==7 in Python (0-based)
47+
# n_res = 2
48+
49+
for jj in range(n_res):
50+
51+
# clear M E - handled by reassignment
52+
53+
if jj == 0:
54+
# M = importdata([outdir,test{k},'_cat_devc.csv'],',',2);
55+
M_data = pd.read_csv(outdir + test[k] + '_cat_devc.csv', skiprows=2)
56+
M = {'data': M_data.values}
57+
58+
# E = importdata([expdir,test{k},'_avg.csv'],',',2);
59+
E_data = pd.read_csv(expdir + test[k] + '_avg.csv', skiprows=2)
60+
E = {'data': E_data.values}
61+
elif jj == 1:
62+
# M = importdata([outdir,test{k},'_fine_cat_devc.csv'],',',2);
63+
M_data = pd.read_csv(outdir + test[k] + '_fine_cat_devc.csv', skiprows=2)
64+
M = {'data': M_data.values}
65+
66+
# E = importdata([expdir,test{k},'_avg.csv'],',',2);
67+
E_data = pd.read_csv(expdir + test[k] + '_avg.csv', skiprows=2)
68+
E = {'data': E_data.values}
69+
70+
# For each experiment, make a contour plot of the extent of a single temperature contour at each time during the experiment
71+
72+
# clear X_mod Y_mod Z_mod X_exp Y_exp Z_exp
73+
74+
# [X_mod,Y_mod] = meshgrid(pos(1:16),M.data(:,1)/60);
75+
X_mod, Y_mod = np.meshgrid(pos[0:16], M['data'][:, 0]/60)
76+
77+
# [X_exp,Y_exp] = meshgrid(pos(1:16),E.data(:,1)/60);
78+
X_exp, Y_exp = np.meshgrid(pos[0:16], E['data'][:, 0]/60)
79+
80+
# Initialize Z_mod and Z_exp
81+
Z_mod = np.zeros((len(M['data'][:, 0]), 16))
82+
Z_exp = np.zeros((len(E['data'][:, 0]), 16))
83+
84+
for kk in range(len(M['data'][:, 0])):
85+
for ii in range(16):
86+
Z_mod[kk, ii] = M['data'][kk, ii+1]
87+
88+
for kk in range(len(E['data'][:, 0])):
89+
for ii in range(16):
90+
Z_exp[kk, ii] = E['data'][kk, ii+4]
91+
92+
newpoints = 100
93+
# [X_mod_interp,Y_mod_interp] = meshgrid(...
94+
# linspace(min(min(X_mod,[],2)),max(max(X_mod,[],2)),newpoints ),...
95+
# linspace(min(min(Y_mod,[],1)),max(max(Y_mod,[],1)),newpoints )...
96+
# );
97+
X_mod_interp, Y_mod_interp = np.meshgrid(
98+
np.linspace(np.min(X_mod), np.max(X_mod), newpoints),
99+
np.linspace(np.min(Y_mod), np.max(Y_mod), newpoints)
100+
)
101+
102+
# Z_mod_interp = interp2(X_mod,Y_mod,Z_mod,X_mod_interp,Y_mod_interp,'makima');
103+
# Flatten the arrays for griddata
104+
points_mod = np.column_stack((X_mod.ravel(), Y_mod.ravel()))
105+
values_mod = Z_mod.ravel()
106+
Z_mod_interp = griddata(points_mod, values_mod,
107+
(X_mod_interp, Y_mod_interp), method='cubic')
108+
109+
# [X_exp_interp,Y_exp_interp] = meshgrid(...
110+
# linspace(min(min(X_exp,[],2)),max(max(X_exp,[],2)),newpoints ),...
111+
# linspace(min(min(Y_exp,[],1)),max(max(Y_exp,[],1)),newpoints )...
112+
# );
113+
X_exp_interp, Y_exp_interp = np.meshgrid(
114+
np.linspace(np.min(X_exp), np.max(X_exp), newpoints),
115+
np.linspace(np.min(Y_exp), np.max(Y_exp), newpoints)
116+
)
117+
118+
# Z_exp_interp = interp2(X_exp,Y_exp,Z_exp,X_exp_interp,Y_exp_interp,'makima');
119+
points_exp = np.column_stack((X_exp.ravel(), Y_exp.ravel()))
120+
values_exp = Z_exp.ravel()
121+
Z_exp_interp = griddata(points_exp, values_exp,
122+
(X_exp_interp, Y_exp_interp), method='cubic')
123+
124+
if jj == 0:
125+
# reset(gca)
126+
# reset(gcf)
127+
plt.figure(figsize=(plot_style["Paper_Width"], plot_style["Paper_Height"]))
128+
mod_symbol = 'r-'
129+
else:
130+
mod_symbol = 'r--'
131+
132+
# [C_mod,h_mod] = contour(X_mod_interp,Y_mod_interp,Z_mod_interp,single_level,mod_symbol) ; hold on
133+
CS_mod = plt.contour(X_mod_interp, Y_mod_interp, Z_mod_interp, single_level, colors='red',
134+
linestyles='-' if jj == 0 else '--')
135+
plt.clabel(CS_mod, fontsize=3, colors='red', inline_spacing=300)
136+
137+
# [C_exp,h_exp] = contour(X_exp_interp,Y_exp_interp,Z_exp_interp,single_level,'k-') ; hold on
138+
CS_exp = plt.contour(X_exp_interp, Y_exp_interp, Z_exp_interp, single_level, colors='black', linestyles='-')
139+
plt.clabel(CS_exp, fontsize=3, colors='black', inline_spacing=300)
140+
141+
# end grid resolution cases
142+
143+
# plot([5.5 5.5],[0 15],'k--')
144+
plt.plot([5.5, 5.5], [0, 15], 'k--')
145+
146+
# plot([0.0 15.],[setpoint(k)/60 setpoint(k)/60],'k:')
147+
plt.plot([0.0, 15.0], [setpoint[k]/60, setpoint[k]/60], 'k:')
148+
149+
# xticks(pos)
150+
# xticklabels({'1.0','1.5','2.0','2.5','3.0','3.5','4.0','4.5','5.0','5.5','6.0','6.5','7.5','9.5','10.5','11.5'})
151+
# ax = gca;
152+
ax = plt.gca()
153+
154+
# ax.XAxis.FontSize = 16;
155+
# ax.YAxis.FontSize = 16;
156+
ax.tick_params(axis='both', which='major', labelsize=16)
157+
158+
# xlabel('Position (m)','FontSize',16,'Interpreter',Font_Interpreter)
159+
plt.xlabel('Position (m)', fontsize=16)
160+
161+
# ylabel('Time (min)','FontSize',16,'Interpreter',Font_Interpreter)
162+
plt.ylabel('Time (min)', fontsize=16)
163+
164+
plt.subplots_adjust(left=plot_style["Plot_X"]/plot_style["Paper_Width"], bottom=plot_style["Plot_Y"]/plot_style["Paper_Height"],
165+
right=(plot_style["Plot_X"]+plot_style["Plot_Width"])/plot_style["Paper_Width"],
166+
top=(plot_style["Plot_Y"]+plot_style["Plot_Height"])/plot_style["Paper_Height"])
167+
168+
plt.rcParams['font.family'] = plot_style["Font_Name"]
169+
170+
plt.axis([0, 15, 0, 15])
171+
172+
plt.text(0.5, 4, test2[k], fontname=plot_style["Font_Name"], fontsize=10)
173+
plt.text(0.5, 3, 'FDS red; Exp black', fontname=plot_style["Font_Name"], fontsize=10)
174+
175+
# Git_Filename = [outdir,test{k},'_cat_git.txt'];
176+
Git_Filename = outdir + test[k] + '_cat_git.txt'
177+
178+
# addverstr(gca,Git_Filename,'linear',0.6,1.05,'Times','TeX',10)
179+
addverstr(ax, Git_Filename, 'linear', 0.6, 1.05, 'Times', 'TeX', 10)
180+
181+
# set(gcf,'Units',Paper_Units);
182+
# set(gcf,'PaperUnits',Paper_Units);
183+
# set(gcf,'PaperSize',[Paper_Width Paper_Height]);
184+
# set(gcf,'Position',[0 0 Paper_Width Paper_Height]);
185+
fig = plt.gcf()
186+
fig.set_size_inches(plot_style["Paper_Width"], plot_style["Paper_Height"])
187+
188+
# print(gcf,'-dpdf',[pltdir,test{k},'_tvT'])
189+
plt.savefig(pltdir + test[k] + '_tvT.pdf', format='pdf', bbox_inches='tight')
190+
191+
# hold off
192+
plt.clf() # Clear the figure for next iteration
193+
194+
# end Experiment loop
195+
196+
print('FHWA_Tunnel completed successfully')
197+
198+

‎Utilities/Python/scripts/UWO_Wind_Tunnel.py‎

Lines changed: 111 additions & 0 deletions
Original file line numberDiff line numberDiff line change
@@ -0,0 +1,111 @@
1+
2+
import os
3+
import csv
4+
import numpy as np
5+
import pandas as pd
6+
7+
# McGrattan
8+
# 07-02-2025
9+
# UWO_Wind Tunnel.py
10+
#
11+
# Reads and transposes line data for UWO Wind Tunnel cases
12+
13+
outdir = '../../../out/UWO_Wind_Tunnel/'
14+
expdir = '../../../exp/UWO_Wind_Tunnel/'
15+
16+
H = [['x', 'Cp_mean', 'Cp_rms', 'Cp_min', 'Cp_max'],
17+
['1.0', 'NaN', 'NaN', 'NaN', 'NaN']]
18+
19+
# Read input CSV file
20+
input_file_path = os.path.join(outdir, 'UWO_inputs.csv')
21+
with open(input_file_path, 'r') as input_file:
22+
csv_reader = csv.reader(input_file)
23+
next(csv_reader) # Skip header line
24+
rows = list(csv_reader)
25+
26+
# Parse the CSV data into columns (equivalent to textscan)
27+
C = [[] for _ in range(38)] # 38 columns based on the textscan format string
28+
29+
for row in rows:
30+
C[0].append(row[0]) # %s - FDS_line_file
31+
C[1].append(row[1]) # %s - output_file
32+
for i in range(2, 38): # remaining columns are integers or floats
33+
if i in [5, 6, 14, 15, 23, 24, 32, 33]: # %f columns
34+
C[i].append(float(row[i]))
35+
else: # %d8 columns
36+
C[i].append(int(row[i]))
37+
38+
FDS_line_file = C[0]
39+
output_file = C[1]
40+
41+
for k in range(len(C[0])):
42+
43+
# Import data (equivalent to importdata)
44+
data_file_path = os.path.join(outdir, C[0][k])
45+
M_data = pd.read_csv(data_file_path, skiprows=2, header=None).values
46+
47+
n_sides = 4
48+
if C[29][k] == 0: # C{30}(k) in MATLAB (0-indexed in Python)
49+
n_sides = 3
50+
if C[20][k] == 0: # C{21}(k) in MATLAB
51+
n_sides = 2
52+
if C[11][k] == 0: # C{12}(k) in MATLAB
53+
n_sides = 1
54+
55+
# Open output file for writing
56+
output_file_path = os.path.join(outdir, output_file[k])
57+
with open(output_file_path, 'w', newline='') as fid:
58+
# Write header
59+
fid.write(','.join(H[0]) + '\n')
60+
61+
dist = 0.0
62+
63+
for j in range(n_sides): # j from 0 to n_sides-1 in Python
64+
65+
coord_col = C[2 + j*9][k] # C{ 3+(j-1)*9}(k) in MATLAB
66+
points = C[3 + j*9][k] # C{ 4+(j-1)*9}(k) in MATLAB
67+
order_index = C[4 + j*9][k] # C{ 5+(j-1)*9}(k) in MATLAB
68+
x1 = C[5 + j*9][k] # C{ 6+(j-1)*9}(k) in MATLAB
69+
x2 = C[6 + j*9][k] # C{ 7+(j-1)*9}(k) in MATLAB
70+
mean_col = C[7 + j*9][k] # C{ 8+(j-1)*9}(k) in MATLAB
71+
rms_col = C[8 + j*9][k] # C{ 9+(j-1)*9}(k) in MATLAB
72+
min_col = C[9 + j*9][k] # C{10+(j-1)*9}(k) in MATLAB
73+
max_col = C[10 + j*9][k] # C{11+(j-1)*9}(k) in MATLAB
74+
75+
# Adjust column indices for 0-based indexing
76+
coord_col -= 1
77+
mean_col -= 1
78+
rms_col -= 1
79+
min_col -= 1
80+
max_col -= 1
81+
82+
if order_index == 1:
83+
for i in range(points): # i from 0 to points-1 in Python
84+
line_data = [
85+
dist + M_data[i, coord_col] - x1,
86+
M_data[i, mean_col],
87+
M_data[i, rms_col],
88+
M_data[i, min_col],
89+
M_data[i, max_col]
90+
]
91+
fid.write('{:6.3f},{:6.3f},{:6.3f},{:6.3f},{:6.3f}\n'.format(*line_data))
92+
else:
93+
for i in range(points-1, -1, -1): # i from points-1 down to 0
94+
line_data = [
95+
dist + x1 - M_data[i, coord_col],
96+
M_data[i, mean_col],
97+
M_data[i, rms_col],
98+
M_data[i, min_col],
99+
M_data[i, max_col]
100+
]
101+
fid.write('{:6.3f},{:6.3f},{:6.3f},{:6.3f},{:6.3f}\n'.format(*line_data))
102+
103+
if j < n_sides - 1: # j+1 < n_sides in MATLAB equivalent
104+
fid.write(','.join(H[1]) + '\n')
105+
106+
dist = dist + abs(x2 - x1)
107+
108+
# Clear M equivalent (Python garbage collection handles this automatically)
109+
del M_data
110+
111+

0 commit comments

Comments
 (0)