import os
import numpy as np
from scipy.integrate import trapezoid
import matplotlib.pyplot as plt
from matplotlib.font_manager import findfont, FontProperties, fontManager
import warnings
from tqdm import tqdm
import pandas as pd
class RadiationCalculator:
"""辐射传热计算类，用于计算样本在不同温度条件下的辐射特性"""
def __init__(self, emissivity_file=r"C:\Users\HP\OneDrive\桌面\样本1发射率.xlsx"):
"""初始化辐射计算器，配置字体和物理常数"""
# 配置中文字体
self._configure_fonts()
# 物理常数
self.h = 6.626e-34 # 普朗克常数 (J·s)
self.c = 2.998e8 # 光速 (m/s)
self.k = 1.381e-23 # 玻尔兹曼常数 (J/K)
# 计算参数
self.lambda_min = 0.3 # 积分下限 (μm)
self.lambda_max = 25.0 # 积分上限 (μm)
# 设置光谱发射率函数
self.set_surface_emissivity_from_excel(emissivity_file)
self.set_atmosphere_emissivity(self.emissivity_a)
# 温度设置
self.object_temp_range = (-50, 100) # 物体温度范围 (°C)
self.object_temp_fixed = 25 # 物体固定温度 (°C)
self.atmosphere_temp_range = (-50, 60) # 大气温度范围 (°C)
self.atmosphere_temp_step = 10 # 大气温度步长 (°C)
self.atmosphere_temp_c = 25 # 大气温度 (°C)，用于物体温度变化模拟
self.temp_step = 10 # 温度步长 (°C)，用于物体温度变化模拟
def _configure_fonts(self):
"""配置matplotlib的中文字体"""
# 备选字体列表，按优先级排序
preferred_fonts = [
'SimHei', # 黑体 - Windows/Linux
'WenQuanYi Micro Hei', # Linux
'Heiti TC', # macOS
'Microsoft YaHei', # Windows
'SimSun', # 宋体 - Windows/Linux
'Heiti TC', # macOS
'Arial Unicode MS', # 支持更多字符
'DejaVu Sans' # 支持更多字符
]
# 尝试指定的备选字体列表
for font in preferred_fonts:
try:
plt.rcParams["font.family"] = font
# 测试字体是否能显示中文
test_text = "测试字体：一阶线性拟合"
fig, ax = plt.subplots(figsize=(1, 1))
ax.text(0.5, 0.5, test_text, fontsize=12, ha='center')
ax.axis('off')
plt.close(fig)
print(f"已成功设置中文字体: {font}")
break
except:
continue
else:
# 如果没有找到中文字体，使用默认字体
plt.rcParams["font.family"] = ["sans-serif"]
print("未找到系统中可用的中文字体，将使用默认字体。图表中的中文可能无法正确显示。")
# 确保负号正确显示
plt.rcParams["axes.unicode_minus"] = False
def blackbody_radiation(self, lambda_m, T):
"""计算黑体辐射强度 (W/(m²·sr·m))"""
lambda_m = np.where(lambda_m <= 0, 1e-10, lambda_m) # 避免波长为零
exponent = self.h * self.c / (lambda_m * self.k * T)
exponent = np.where(exponent > 709, 709, exponent) # 避免数值溢出
return (2 * self.h * self.c**2 / lambda_m**5) / (np.exp(exponent) - 1)
# 大气光谱发射率函数
def emissivity_a(self, lambda_m):
"""大气的光谱发射率函数"""
lambda_um = lambda_m * 1e6
emissivity = np.piecewise(lambda_um,
[lambda_um < 0, (lambda_um >= 0) & (lambda_um < 8),
(lambda_um >= 8) & (lambda_um <= 14), lambda_um > 14],
[0, 1.0, 0.52, 1.0])
return emissivity
def set_surface_emissivity_from_excel(self, file_path):
"""从Excel文件中设置物体表面的光谱发射率函数"""
data = pd.read_excel(file_path)
wavelengths = data['波长'].values * 1e-6 # 转换为米
emissivities = data['发射率'].values
# 创建插值函数
self.surface_emissivity_func = lambda x: np.interp(x, wavelengths, emissivities, left=0, right=0)
def set_atmosphere_emissivity(self, emissivity_func):
"""设置大气的光谱发射率函数"""
self.atmosphere_emissivity_func = emissivity_func
def calculate_radiation_power(self, T_s, T_a):
"""计算物体的辐射功率、吸收的大气辐射功率和净辐射功率"""
lambda_min_m = self.lambda_min * 1e-6
lambda_max_m = self.lambda_max * 1e-6
# 获取实际的波长范围内的数据点
data = pd.read_excel(r"C:\Users\HP\OneDrive\桌面\样本1发射率.xlsx")
wavelengths_m = data['波长'].values * 1e-6 # 转换为米
emissivities = data['发射率'].values
# 过滤掉不在积分范围内的数据点
mask = (wavelengths_m >= lambda_min_m) & (wavelengths_m <= lambda_max_m)
wavelengths_m = wavelengths_m[mask]
emissivities = emissivities[mask]
if len(wavelengths_m) == 0:
return 0, 0, 0
# 计算物体辐射功率
radiance_s = self.blackbody_radiation(wavelengths_m, T_s)
P_rad = np.pi * trapezoid(radiance_s * emissivities, wavelengths_m)
# 计算吸收的大气辐射功率
radiance_a = self.blackbody_radiation(wavelengths_m, T_a)
emissivity_a = self.atmosphere_emissivity_func(wavelengths_m)
P_atmosphere = np.pi * trapezoid(radiance_a * emissivity_a * emissivities, wavelengths_m)
# 计算净辐射功率
P_net = P_rad - P_atmosphere
return P_rad, P_atmosphere, P_net
def run_simulation_object_temperature(self):
"""运行模拟过程，计算样本在不同物体温度下的辐射特性"""
object_temp_min_k = self.object_temp_range[0] + 273.15
object_temp_max_k = self.object_temp_range[1] + 273.15
object_temps_k = np.arange(object_temp_min_k, object_temp_max_k + 1, self.temp_step)
object_temps_c = object_temps_k - 273.15
atmosphere_temp_k = self.atmosphere_temp_c + 273.15
print(f"物体温度范围: {min(object_temps_c):.1f}°C 到 {max(object_temps_c):.1f}°C")
print(f"大气温度固定为: {self.atmosphere_temp_c}°C")
print(f"温度步长: {self.temp_step}°C")
results = {
'object_temp_c': [],
'net_power': []
}
total_iterations = len(object_temps_k)
progress_bar = tqdm(total=total_iterations, desc="计算辐射功率")
for T_s_k in object_temps_k:
T_s_c = T_s_k - 273.15
_, _, P_net = self.calculate_radiation_power(T_s_k, atmosphere_temp_k)
results['object_temp_c'].append(float(T_s_c))
results['net_power'].append(float(P_net))
progress_bar.update(1)
progress_bar.close()
return pd.DataFrame(results)
def run_simulation_atmosphere_temperature(self):
"""运行模拟过程，计算样本在不同大气温度下的辐射特性"""
atmosphere_temp_min_k = self.atmosphere_temp_range[0] + 273.15
atmosphere_temp_max_k = self.atmosphere_temp_range[1] + 273.15
atmosphere_temps_k = np.arange(atmosphere_temp_min_k, atmosphere_temp_max_k + 1, self.atmosphere_temp_step)
atmosphere_temps_c = atmosphere_temps_k - 273.15
object_temp_k = self.object_temp_fixed + 273.15
print(f"大气温度范围: {min(atmosphere_temps_c):.1f}°C 到 {max(atmosphere_temps_c):.1f}°C")
print(f"物体温度固定为: {self.object_temp_fixed}°C")
print(f"温度步长: {self.atmosphere_temp_step}°C")
results = {
'atmosphere_temp_c': [],
'net_power': []
}
total_iterations = len(atmosphere_temps_k)
progress_bar = tqdm(total=total_iterations, desc="计算辐射功率")
for T_a_k in atmosphere_temps_k:
T_a_c = T_a_k - 273.15
_, _, P_net = self.calculate_radiation_power(object_temp_k, T_a_k)
results['atmosphere_temp_c'].append(float(T_a_c))
results['net_power'].append(float(P_net))
progress_bar.update(1)
progress_bar.close()
return pd.DataFrame(results)
def run_simulation_double_temperature(self):
"""运行模拟过程，计算样本在不同物体温度和大气温度组合下的辐射特性"""
object_temp_min_k = self.object_temp_range[0] + 273.15
object_temp_max_k = self.object_temp_range[1] + 273.15
object_temps_k = np.arange(object_temp_min_k, object_temp_max_k + 1, self.temp_step)
object_temps_c = object_temps_k - 273.15
atmosphere_temp_min_k = self.atmosphere_temp_range[0] + 273.15
atmosphere_temp_max_k = self.atmosphere_temp_range[1] + 273.15
atmosphere_temps_k = np.arange(atmosphere_temp_min_k, atmosphere_temp_max_k + 1, self.atmosphere_temp_step)
atmosphere_temps_c = atmosphere_temps_k - 273.15
print(f"物体温度范围: {min(object_temps_c):.1f}°C 到 {max(object_temps_c):.1f}°C")
print(f"大气温度范围: {min(atmosphere_temps_c):.1f}°C 到 {max(atmosphere_temps_c):.1f}°C")
print(f"温度步长: {self.temp_step}°C 和 {self.atmosphere_temp_step}°C")
results = {
'object_temp_c': [],
'atmosphere_temp_c': [],
'net_power': []
}
total_iterations = len(object_temps_k) * len(atmosphere_temps_k)
progress_bar = tqdm(total=total_iterations, desc="计算辐射功率")
for T_s_k in object_temps_k:
for T_a_k in atmosphere_temps_k:
T_s_c = T_s_k - 273.15
T_a_c = T_a_k - 273.15
_, _, P_net = self.calculate_radiation_power(T_s_k, T_a_k)
results['object_temp_c'].append(float(T_s_c))
results['atmosphere_temp_c'].append(float(T_a_c))
results['net_power'].append(float(P_net))
progress_bar.update(1)
progress_bar.close()
return pd.DataFrame(results)
def save_results_to_excel(self, results_obj_temp, results_atmo_temp, results_double_temp, file_path):
"""将结果保存到Excel文件的不同工作表中"""
# 确保输出目录存在
output_dir = os.path.dirname(file_path)
if not os.path.exists(output_dir):
os.makedirs(output_dir)
print(f"创建输出目录: {output_dir}")
with pd.ExcelWriter(file_path) as writer:
results_obj_temp.to_excel(writer, sheet_name='Object_Temperature_Change', index=False)
results_atmo_temp.to_excel(writer, sheet_name='Atmosphere_Temperature_Change', index=False)
results_double_temp.to_excel(writer, sheet_name='Double_Temperature_Change', index=False)
print(f"结果已保存到 {file_path}")
def plot_and_save_figures(self, results_obj_temp, results_atmo_temp, results_double_temp, output_folder):
"""绘制并保存图表"""
# 确保输出目录存在
if not os.path.exists(output_folder):
os.makedirs(output_folder)
print(f"创建输出目录: {output_folder}")
# 绘制物体温度变化图
plt.figure(figsize=(10, 6))
plt.scatter(results_obj_temp['object_temp_c'], results_obj_temp['net_power'], label='数据点')
# 一次函数拟合
coeffs_obj_temp = np.polyfit(results_obj_temp['object_temp_c'], results_obj_temp['net_power'], 1)
poly_fit_obj_temp = np.poly1d(coeffs_obj_temp)
plt.plot(results_obj_temp['object_temp_c'], poly_fit_obj_temp(results_obj_temp['object_temp_c']), color='red', label=f'y={coeffs_obj_temp[0]:.2f}x+{coeffs_obj_temp[1]:.2f}')
# 计算R²值
residuals_obj_temp = results_obj_temp['net_power'] - poly_fit_obj_temp(results_obj_temp['object_temp_c'])
ss_res_obj_temp = np.sum(residuals_obj_temp**2)
ss_tot_obj_temp = np.sum((results_obj_temp['net_power'] - np.mean(results_obj_temp['net_power']))**2)
r_squared_obj_temp = 1 - (ss_res_obj_temp / ss_tot_obj_temp)
# 计算绝对误差和相对误差
abs_errors_obj_temp = np.abs(results_obj_temp['net_power'] - poly_fit_obj_temp(results_obj_temp['object_temp_c']))
rel_errors_obj_temp = abs_errors_obj_temp / np.abs(results_obj_temp['net_power']) * 100 # 相对误差转为百分比
# 计算平均绝对误差和平均相对误差
mae_obj_temp = np.mean(abs_errors_obj_temp)
mre_obj_temp = np.mean(rel_errors_obj_temp)
# 计算均方根相对误差 (RMSRE)
rmsre_obj_temp = np.sqrt(np.mean((rel_errors_obj_temp / 100) ** 2)) * 100
plt.xlabel('物体温度 (°C)')
plt.ylabel('净辐射功率 (W/m²)')
plt.title(f'在固定大气温度为25℃下不同物体温度对净辐射功率的影响\nR²={r_squared_obj_temp:.2f}, MAE={mae_obj_temp:.2f} W/m², MRE={mre_obj_temp:.2f}%, RMSRE={rmsre_obj_temp:.2f}%')
plt.legend(loc='upper right')
plt.grid(True)
obj_temp_plot_path = os.path.join(output_folder, 'object_temperature_plot.svg')
plt.savefig(obj_temp_plot_path, format='svg')
plt.show()
# 绘制物体温度变化图的绝对误差
plt.figure(figsize=(10, 6))
plt.scatter(results_obj_temp['object_temp_c'], abs_errors_obj_temp, label='绝对误差')
plt.axhline(y=mae_obj_temp, color='green', linestyle='--', label=f'平均绝对误差 ({mae_obj_temp:.2f} W/m²)')
plt.xlabel('物体温度 (°C)')
plt.ylabel('绝对误差 (W/m²)')
plt.title(f'在固定大气温度为25℃下不同物体温度对净辐射功率的绝对误差')
plt.legend(loc='upper right')
plt.grid(True)
obj_abs_error_plot_path = os.path.join(output_folder, 'object_absolute_error_plot.svg')
plt.savefig(obj_abs_error_plot_path, format='svg')
plt.show()
# 绘制物体温度变化图的相对误差
plt.figure(figsize=(10, 6))
plt.scatter(results_obj_temp['object_temp_c'], rel_errors_obj_temp, label='相对误差 (%)')
plt.axhline(y=mre_obj_temp, color='red', linestyle='--', label=f'平均相对误差 ({mre_obj_temp:.2f}%)')
plt.axhline(y=rmsre_obj_temp, color='purple', linestyle='-.', label=f'RMSRE ({rmsre_obj_temp:.2f}%)')
plt.xlabel('物体温度 (°C)')
plt.ylabel('相对误差 (%)')
plt.title(f'在固定大气温度为25℃下不同物体温度对净辐射功率的相对误差')
plt.legend(loc='upper right')
plt.grid(True)
obj_rel_error_plot_path = os.path.join(output_folder, 'object_relative_error_plot.svg')
plt.savefig(obj_rel_error_plot_path, format='svg')
plt.show()
# 绘制大气温度变化图
plt.figure(figsize=(10, 6))
plt.scatter(results_atmo_temp['atmosphere_temp_c'], results_atmo_temp['net_power'], label='数据点')
# 一次函数拟合
coeffs_atmo_temp = np.polyfit(results_atmo_temp['atmosphere_temp_c'], results_atmo_temp['net_power'], 1)
poly_fit_atmo_temp = np.poly1d(coeffs_atmo_temp)
plt.plot(results_atmo_temp['atmosphere_temp_c'], poly_fit_atmo_temp(results_atmo_temp['atmosphere_temp_c']), color='red', label=f'y={coeffs_atmo_temp[0]:.2f}x+{coeffs_atmo_temp[1]:.2f}')
# 计算R²值
residuals_atmo_temp = results_atmo_temp['net_power'] - poly_fit_atmo_temp(results_atmo_temp['atmosphere_temp_c'])
ss_res_atmo_temp = np.sum(residuals_atmo_temp**2)
ss_tot_atmo_temp = np.sum((results_atmo_temp['net_power'] - np.mean(results_atmo_temp['net_power']))**2)
r_squared_atmo_temp = 1 - (ss_res_atmo_temp / ss_tot_atmo_temp)
# 计算绝对误差和相对误差
abs_errors_atmo_temp = np.abs(results_atmo_temp['net_power'] - poly_fit_atmo_temp(results_atmo_temp['atmosphere_temp_c']))
rel_errors_atmo_temp = abs_errors_atmo_temp / np.abs(results_atmo_temp['net_power']) * 100 # 相对误差转为百分比
# 计算平均绝对误差和平均相对误差
mae_atmo_temp = np.mean(abs_errors_atmo_temp)
mre_atmo_temp = np.mean(rel_errors_atmo_temp)
# 计算均方根相对误差 (RMSRE)
rmsre_atmo_temp = np.sqrt(np.mean((rel_errors_atmo_temp / 100) ** 2)) * 100
plt.xlabel('大气温度 (°C)')
plt.ylabel('净辐射功率 (W/m²)')
plt.title(f'在固定物体温度为25℃下不同大气温度对净辐射功率的影响\nR²={r_squared_atmo_temp:.2f}, MAE={mae_atmo_temp:.2f} W/m², MRE={mre_atmo_temp:.2f}%, RMSRE={rmsre_atmo_temp:.2f}%')
plt.legend(loc='upper right')
plt.grid(True)
atmo_temp_plot_path = os.path.join(output_folder, 'atmosphere_temperature_plot.svg')
plt.savefig(atmo_temp_plot_path, format='svg')
plt.show()
# 绘制大气温度变化图的绝对误差
plt.figure(figsize=(10, 6))
plt.scatter(results_atmo_temp['atmosphere_temp_c'], abs_errors_atmo_temp, label='绝对误差')
plt.axhline(y=mae_atmo_temp, color='green', linestyle='--', label=f'平均绝对误差 ({mae_atmo_temp:.2f} W/m²)')
plt.xlabel('大气温度 (°C)')
plt.ylabel('绝对误差 (W/m²)')
plt.title(f'在固定物体温度为25℃下不同大气温度对净辐射功率的绝对误差')
plt.legend(loc='upper right')
plt.grid(True)
atmo_abs_error_plot_path = os.path.join(output_folder, 'atmosphere_absolute_error_plot.svg')
plt.savefig(atmo_abs_error_plot_path, format='svg')
plt.show()
# 绘制大气温度变化图的相对误差
plt.figure(figsize=(10, 6))
plt.scatter(results_atmo_temp['atmosphere_temp_c'], rel_errors_atmo_temp, label='相对误差 (%)')
plt.axhline(y=mre_atmo_temp, color='red', linestyle='--', label=f'平均相对误差 ({mre_atmo_temp:.2f}%)')
plt.axhline(y=rmsre_atmo_temp, color='purple', linestyle='-.', label=f'RMSRE ({rmsre_atmo_temp:.2f}%)')
plt.xlabel('大气温度 (°C)')
plt.ylabel('相对误差 (%)')
plt.title(f'在固定物体温度为25℃下不同大气温度对净辐射功率的相对误差')
plt.legend(loc='upper right')
plt.grid(True)
atmo_rel_error_plot_path = os.path.join(output_folder, 'atmosphere_relative_error_plot.svg')
plt.savefig(atmo_rel_error_plot_path, format='svg')
plt.show()
# 绘制双温度变量的变化图
plt.figure(figsize=(10, 6))
sc = plt.scatter(results_double_temp['atmosphere_temp_c'], results_double_temp['object_temp_c'], c=results_double_temp['net_power'], cmap='viridis', s=50)
plt.colorbar(sc, label='净辐射功率 (W/m²)', pad=0.01)
plt.xlabel('大气温度 (°C)')
plt.ylabel('物体温度 (°C)')
plt.title('不同物体温度和大气温度组合下的净辐射功率')
plt.grid(True)
double_temp_scatter_path = os.path.join(output_folder, 'double_temperature_scatter_plot.svg')
plt.savefig(double_temp_scatter_path, format='svg')
plt.show()
# 线性拟合
X = np.column_stack((results_double_temp['atmosphere_temp_c'], results_double_temp['object_temp_c']))
coeffs_double_temp, residuals, rank, s = np.linalg.lstsq(X, results_double_temp['net_power'], rcond=None)
# 打印系数以调试
print(f"拟合得到的系数: {coeffs_double_temp}")
print(f"残差平方和: {residuals}")
print(f"秩: {rank}")
print(f"singular values: {s}")
# 确保 coeffs_double_temp 是一个 numpy 数组
coeffs_double_temp = np.array(coeffs_double_temp)
if len(coeffs_double_temp) == 2:
a, b = coeffs_double_temp
c = 0 # 假设截距为0
elif len(coeffs_double_temp) == 3:
a, b = coeffs_double_temp[:2]
c = coeffs_double_temp[2]
else:
raise ValueError("拟合得到的系数数量不符合预期")
# 计算预测值
predicted_values = a * results_double_temp['atmosphere_temp_c'].to_numpy() + b * results_double_temp['object_temp_c'].to_numpy() + c
# 计算绝对误差和相对误差
abs_errors_double_temp = np.abs(results_double_temp['net_power'] - predicted_values)
rel_errors_double_temp = abs_errors_double_temp / np.abs(results_double_temp['net_power']) * 100 # 相对误差转为百分比
# 计算均方根相对误差 (RMSRE)
rmsre_double_temp = np.sqrt(np.mean((rel_errors_double_temp / 100) ** 2)) * 100
# 筛选相对误差大于500%的点
mask_high_error = rel_errors_double_temp > 500
high_error_points = results_double_temp[mask_high_error]
high_error_rel_errors = rel_errors_double_temp[mask_high_error]
# 绘制绝对误差图
plt.figure(figsize=(10, 6))
# 正常点
plt.scatter(
results_double_temp['atmosphere_temp_c'][~mask_high_error],
results_double_temp['object_temp_c'][~mask_high_error],
c=abs_errors_double_temp[~mask_high_error],
cmap='coolwarm',
s=50,
label='相对误差 ≤ 500%'
)
# 特殊点（相对误差 > 500%）
plt.scatter(
high_error_points['atmosphere_temp_c'],
high_error_points['object_temp_c'],
c='black',
s=100,
marker='x',
label=f'相对误差 > 500% ({len(high_error_points)}个点)'
)
plt.colorbar(label='绝对误差 (W/m²)', pad=0.01)
plt.xlabel('大气温度 (°C)')
plt.ylabel('物体温度 (°C)')
plt.title(f'不同物体温度和大气温度组合下的绝对误差\nRMSRE={rmsre_double_temp:.2f}%')
plt.legend(loc='upper right')
plt.grid(True)
double_temp_abs_error_path = os.path.join(output_folder, 'double_temperature_absolute_error_plot.svg')
plt.savefig(double_temp_abs_error_path, format='svg')
plt.show()
# 绘制相对误差图（使用更鲜明的红色调，限制颜色映射范围为0-500%）
plt.figure(figsize=(10, 6))
# 创建更鲜明的红色调的colormap（使用_r反转以获得更暗的红色）
red_cmap = plt.cm.get_cmap('Reds_r')
# 正常点（限制颜色映射范围为0-500%，增加点大小）
scatter = plt.scatter(
results_double_temp['atmosphere_temp_c'][~mask_high_error],
results_double_temp['object_temp_c'][~mask_high_error],
c=np.clip(rel_errors_double_temp[~mask_high_error], 0, 500),
cmap=red_cmap,
s=60, # 增加点大小，提高可见性
alpha=0.8 # 增加不透明度，使颜色更鲜明
)
# 特殊点（相对误差 > 500%，增加点大小）
plt.scatter(
high_error_points['atmosphere_temp_c'],
high_error_points['object_temp_c'],
c='black',
s=120, # 增加特殊点大小，使其更突出
marker='x',
linewidth=2, # 增加线条宽度
label=f'相对误差 > 500% ({len(high_error_points)}个点)'
)
# 添加颜色条，并设置为百分比格式
cbar = plt.colorbar(scatter, pad=0.01)
cbar.set_label('相对误差 (%)', rotation=270, labelpad=20)
# 设置标题和标签
plt.xlabel('大气温度 (°C)')
plt.ylabel('物体温度 (°C)')
plt.title(f'不同物体温度和大气温度组合下的相对误差\nRMSRE={rmsre_double_temp:.2f}%')
plt.legend(loc='upper right')
plt.grid(True, linestyle='--', alpha=0.6) # 优化网格线
# 保存并显示图表
double_temp_rel_error_path = os.path.join(output_folder, 'double_temperature_relative_error_plot.svg')
plt.savefig(double_temp_rel_error_path, format='svg')
plt.show()
# 绘制拟合平面
xx, yy = np.meshgrid(np.linspace(min(results_double_temp['atmosphere_temp_c']), max(results_double_temp['atmosphere_temp_c']), 10),
np.linspace(min(results_double_temp['object_temp_c']), max(results_double_temp['object_temp_c']), 10))
zz = a * xx + b * yy + c
plt.figure(figsize=(10, 6))
ax = plt.axes(projection='3d')
# 正常点
scatter = ax.scatter(
results_double_temp['atmosphere_temp_c'][~mask_high_error],
results_double_temp['object_temp_c'][~mask_high_error],
results_double_temp['net_power'][~mask_high_error],
c=results_double_temp['net_power'][~mask_high_error],
cmap='viridis',
s=50,
label='相对误差 ≤ 500%'
)
# 特殊点（相对误差 > 500%）
ax.scatter(
high_error_points['atmosphere_temp_c'],
high_error_points['object_temp_c'],
high_error_points['net_power'],
c='black',
s=100,
marker='x',
label=f'相对误差 > 500% ({len(high_error_points)}个点)'
)
surface = ax.plot_surface(xx, yy, zz, alpha=0.5, cmap='coolwarm')
ax.set_xlabel('大气温度 (°C)')
ax.set_ylabel('物体温度 (°C)')
ax.set_zlabel('净辐射功率 (W/m²)')
ax.zaxis.labelpad = 20
ax.xaxis.labelpad = 20
ax.yaxis.labelpad = 20
ax.set_title(f'不同物体温度和大气温度组合下的净辐射功率\n拟合方程: z={a:.2f}x + {b:.2f}y + {c:.2f}\nRMSRE={rmsre_double_temp:.2f}%')
plt.colorbar(scatter, label='净辐射功率 (W/m²)', ax=ax, location='left')
plt.colorbar(surface, label='拟合面', ax=ax, location='right')
plt.legend(loc='upper right')
plt.grid(True)
double_temp_fit_path = os.path.join(output_folder, 'double_temperature_fit_plot.svg')
plt.savefig(double_temp_fit_path, format='svg')
plt.show()
# 示例用法
if __name__ == "__main__":
calculator = RadiationCalculator()
# 运行模拟，计算物体温度变化的结果
results_obj_temp = calculator.run_simulation_object_temperature()
# 运行模拟，计算大气温度变化的结果
results_atmo_temp = calculator.run_simulation_atmosphere_temperature()
# 运行模拟，计算双温度变量变化的结果
results_double_temp = calculator.run_simulation_double_temperature()
# 创建输出文件夹
output_folder = r"D:/代码/梯形法则实际状况"
# 保存结果到Excel文件
excel_file_path = os.path.join(output_folder, 'radiation_results.xlsx')
calculator.save_results_to_excel(results_obj_temp, results_atmo_temp, results_double_temp, excel_file_path)
# 绘制并保存图表
calculator.plot_and_save_figures(results_obj_temp, results_atmo_temp, results_double_temp, output_folder)这个代码的计算模块计算是否出现问题了
扫一扫
9404
被折叠的 条评论
为什么被折叠?
到【灌水乐园】发言
