主题027：全谱线-by-线计算

摘要

全谱线-by-线（Line-by-Line, LBL）计算是气体辐射传热中精度最高的计算方法。本教程系统介绍LBL计算的理论基础、实施方法和工程应用。通过详细讲解HITRAN和HITEMP等光谱数据库的使用方法，以及高分辨率光谱计算的关键技术，帮助读者掌握LBL计算的核心原理。结合丰富的Python仿真案例，展示如何在实际工程问题中应用LBL计算作为基准解，验证其他简化模型的精度。

关键词

线-by-线计算、高分辨率光谱、HITRAN、HITEMP、光谱数据库、逐线积分、Voigt线型、辐射特性

---

1. 引言

1.1 LBL计算的定位

在气体辐射计算方法的精度谱系中，LBL计算处于最高精度端：

精度层次：

1. 逐线计算（LBL）：最高精度，计算每条谱线的详细贡献
2. 窄谱带模型（NBM）：中等精度，统计描述谱线群
3. 宽谱带模型（WBM）：较低精度，整体描述振动-转动带
4. 全谱模型（WSGGM等）：最低精度，全光谱平均

计算成本对比（相对值）：

• LBL：1000-10000
• NBM：50-100
• WBM：5-10
• WSGGM：1

1.2 LBL计算的重要性

尽管计算成本高昂，LBL计算在以下场景不可或缺：

1. 基准验证

• 作为"数值实验"验证简化模型
• 评估谱带模型的适用范围
• 确定模型参数

2. 高精度需求场景

• 大气遥感反演
• 高精度光谱测量
• 科学研究

3. 极端条件

• 高温（>2000K）
• 高压（>10 atm）
• 非平衡流动

1.3 LBL计算的发展历史

早期阶段（1960s-1980s）：

• 基于实验室测量的线参数表
• 有限的谱线数据
• 计算能力受限

数据库时代（1990s-2000s）：

• HITRAN数据库建立（Rothman et al., 1992）
• HITEMP高温数据库发布
• 计算能力大幅提升

现代发展（2010s至今）：

• 数据库持续更新（HITRAN2020）
• GPU加速计算
• 与CFD的深度集成

---

2. 光谱数据库

2.1 HITRAN数据库

**HITRAN（High-resolution TRANsmission）**是目前最广泛使用的分子光谱数据库。

基本信息：

• 维护机构：哈佛-史密松天体物理中心（CfA）
• 首发时间：1973年
• 最新版本：HITRAN2020
• 网址：https://hitran.org

数据内容：

• 49种分子
• 超过500万条谱线
• 温度范围：主要适用于T < 1000K
• 波数范围：0-35000 cm⁻¹

数据格式：
每条谱线包含以下参数：

ν₀: 谱线中心波数 [cm⁻¹]
S: 谱线强度 [cm⁻¹/(molecule·cm⁻²)]
γ_air: 空气展宽半宽度 [cm⁻¹/atm]
γ_self: 自展宽半宽度 [cm⁻¹/atm]
E_low: 低能级能量 [cm⁻¹]
n: 温度依赖指数
δ: 压力位移 [cm⁻¹/atm]

2.2 HITEMP数据库

HITEMP是专门为高温应用开发的数据库。

基本信息：

• 维护机构：与HITRAN相同
• 首发时间：2010年
• 最新版本：HITEMP2010

与HITRAN的区别：

• 包含更多高温激发的谱线
• 适用于T > 1000K的条件
• 目前仅包含H₂O、CO₂、CO、NO等主要分子
• 数据量更大（H₂O超过1亿条谱线）

温度适用范围：

• H₂O: up to 3000K
• CO₂: up to 2000K
• CO: up to 5000K

2.3 其他数据库

GEISA（Gestion et Etude des Informations Spectroscopiques Atmosphériques）：

• 法国开发
• 主要用于大气研究
• 与HITRAN类似

CDSD（Carbon Dioxide Spectroscopic Databank）：

• 专门针对CO₂
• 高温精度更高
• 适用于燃烧应用

NIST化学数据库：

• 包含原子谱线数据
• 用于等离子体辐射

2.4 数据库的使用方法

数据下载：

# 使用HAPI (HITRAN Application Programming Interface)
from hapi import *

# 获取数据
downloadMolecule('H2O', isotopologue='all')
downloadMolecule('CO2', isotopologue='all')

# 读取数据
fetch('H2O_1200K', 1, 1, 1000, 3000)  # H2O, 1000-3000 cm⁻¹

数据查询：

# 查询特定波数范围的谱线
nu, sw = getColumns('H2O_1200K', ['nu', 'sw'])
mask = (nu > 2000) & (nu < 2500)
nu_selected = nu[mask]
sw_selected = sw[mask]

---

3. LBL计算的理论基础

3.1 单条谱线的吸收系数

洛伦兹线型（碰撞展宽主导）：

$${\kappa }_L(\nu )=\frac{S}{\pi }\frac{{\gamma }_L}{(\nu -{\nu }_0{)}^2+{\gamma }_L^2}$$

其中：

• $S$: 谱线强度
• ${\gamma }_L$: 洛伦兹半宽度
• ${\nu }_0$: 谱线中心波数

多普勒线型（热运动展宽主导）：

$${\kappa }_D(\nu )=\frac{S}{{\gamma }_D}\sqrt{\frac{\ln 2}{\pi }}\exp \left(-\ln 2{\left(\frac{\nu -{\nu }_0}{{\gamma }_D}\right)}^2\right)$$

其中多普勒半宽度：

$${\gamma }_D=\frac{{\nu }_0}{c}\sqrt{\frac{2kT\ln 2}{m}}$$

Voigt线型（两种展宽都重要）：

Voigt线型是洛伦兹线型和多普勒线型的卷积：

$${\kappa }_V(\nu )={\int }_{-\infty }^{\infty }{\kappa }_L({\nu }^{\prime }){\kappa }_D(\nu -{\nu }^{\prime })d{\nu }^{\prime }$$

实际计算使用近似公式或数值积分。

3.2 谱线强度的温度依赖

谱线强度随温度变化：

$$S(T)=S(T_{ref})\frac{Q(T_{ref})}{Q(T)}\frac{\exp (-c_2{E}_{low}/T)}{\exp (-c_2{E}_{low}/T_{ref})}\frac{1-\exp (-c_2{\nu }_0/T)}{1-\exp (-c_2{\nu }_0/T_{ref})}$$

其中：

• $Q(T)$: 配分函数
• $E_{low}$: 低能级能量
• $c_2=1.4388$ cm·K（第二辐射常数）

3.3 半宽度的温度和压力依赖

洛伦兹半宽度：

$${\gamma }_L={\left(\frac{T_{ref}}{T}\right)}^n\left[{\gamma }_{air}(P-P_s)+{\gamma }_{self}{P}_s\right]$$

其中：

• $n$: 温度依赖指数（通常0.5-0.75）
• $P_s$: 分压
• ${\gamma }_{air}$, ${\gamma }_{self}$: 空气和自展宽系数

多普勒半宽度：

$${\gamma }_D\propto \sqrt{T}$$

3.4 光谱吸收系数的计算

在波数$\nu$处的总吸收系数是所有谱线贡献的叠加：

$$\kappa (\nu )=\sum _i{\kappa }_i(\nu )$$

计算策略：

1. 直接求和：对所有谱线直接计算并叠加

   • 精度最高
   • 计算量最大

2. 截断处理：只考虑中心在一定范围内的谱线

   • 设置截断距离（通常25-100 cm⁻¹）
   • 平衡精度和效率

3. 线翼截断：对远翼使用近似

   • 远翼贡献较小
   • 可使用子洛伦兹线型

---

4. LBL计算的实现

4.1 计算流程

开始
  ↓
读取光谱数据库
  ↓
设置计算参数（T, P, 组分）
  ↓
选择波数网格
  ↓
对每个波数点：
  ├── 计算每条谱线的贡献
  ├── 叠加所有贡献
  └── 得到κ(ν)
  ↓
计算辐射特性（发射率、透射率）
  ↓
输出结果

4.2 波数网格的选择

均匀网格：

• 最简单
• 需要非常细的网格（~0.001 cm⁻¹）
• 计算量大

自适应网格：

• 在谱线密集区加密
• 在谱线稀疏区稀疏
• 提高效率

谱线中心网格：

• 只在谱线中心附近加密
• 线翼使用解析近似
• 工程常用

4.3 计算优化技术

1. 谱线分组：

• 将相近的谱线分组处理
• 减少重复计算

2. 快速Voigt计算：

• 使用近似公式（如Humlicek算法）
• 查表法
• 避免数值积分

3. 并行计算：

• 波数并行：不同波数点独立计算
• 谱线并行：不同谱线独立计算
• GPU加速

4. 预计算：

• 预先计算线型函数
• 建立查找表

4.4 与RTE求解的耦合

直接耦合：

• 每个RTE求解步骤都进行LBL计算
• 精度最高，计算成本不可接受

查表法：

• 预先在(T, P, L)空间计算
• RTE求解时插值
• 工程实用

混合方法：

• 关键区域使用LBL
• 其他区域使用简化模型

---

5. 高分辨率计算的关键技术

5.1 线翼修正

标准洛伦兹线型在远翼（|ν-ν₀| >> γ）可能过高估计吸收。

子洛伦兹线型（Sub-Lorentzian）：

$${\kappa }_{far-wing}={\kappa }_L\cdot \chi (|\nu -{\nu }_0|)$$

其中χ是修正函数，通常<1。

连续吸收：

• 由远翼叠加形成
• 在窗口区域重要
• HITRAN提供连续吸收系数

5.2 谱线混合

在高压条件下，谱线可能相互影响（混合效应）。

Rosenkranz近似：

$${\kappa }_{mixed}={\kappa }_{isolated}+\Delta \kappa$$

其中Δκ是混合修正项。

5.3 非洛伦兹线型

某些情况下需要使用更复杂的线型：

Galatry线型：考虑分子间碰撞
Rautian线型：考虑速度变化
Speed-dependent Voigt：考虑速度依赖的展宽

5.4 分子间相互作用

CO₂-H₂O相互作用：

• 形成复合物
• 产生额外的连续吸收
• 在高温高压下重要

压力诱导吸收：

• 由分子间碰撞引起
• 在红外窗口区域显著
• 需要特殊处理

---

6. LBL计算的精度与验证

6.1 误差来源分析

1. 数据库误差：

• 谱线参数测量误差
• 缺失谱线（特别是弱线）
• 高温外推误差

2. 计算误差：

• 波数网格不够密
• 线翼截断
• 数值舍入

3. 物理模型误差：

• 线型选择
• 线混合处理
• 连续吸收

6.2 与实验对比

实验室测量：

• 傅里叶变换红外光谱（FTIR）
• 可调谐激光吸收光谱（TDLAS）
• 腔衰荡光谱（CRDS）

对比指标：

• 透射率误差 < 1%
• 发射率误差 < 2%
• 谱线位置误差 < 0.001 cm⁻¹

6.3 不确定性量化

蒙特卡洛方法：

• 对输入参数进行随机扰动
• 统计输出结果分布
• 评估不确定性

敏感性分析：

• 识别关键参数
• 指导数据库改进
• 优化计算策略

---

7. 工程应用

7.1 大气辐射传输

应用场景：

• 卫星遥感反演
• 气候模式
• 大气校正

特点：

• 路径长（几十公里）
• 分层计算
• 需要多次散射

常用软件：

• LBLRTM（Line-By-Line Radiative Transfer Model）
• GENLN2
• ARTS

7.2 燃烧诊断

应用场景：

• 温度测量（TDLAS）
• 浓度测量
• 流速测量

特点：

• 高温（1000-3000K）
• 使用HITEMP数据库
• 需要快速计算

7.3 基准计算

应用场景：

• 验证WSGGM参数
• 评估谱带模型精度
• CFD辐射模型验证

典型设置：

• 一维等温层
• 均匀或非均匀路径
• 对比发射率、热流

---

8. Python仿真案例

以下是基于LBL计算的Python仿真程序，包含多个实际应用案例。

案例1：单条谱线线型分析

"""
案例1：单条谱线线型分析
对比洛伦兹、多普勒和Voigt线型
"""
import numpy as np
import matplotlib.pyplot as plt
from scipy.special import wofz
import matplotlib
matplotlib.use('Agg')

# 设置中文字体
plt.rcParams['font.sans-serif'] = ['SimHei', 'DejaVu Sans', 'Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False

def lorentz_line(nu, nu0, S, gamma_L):
    """洛伦兹线型"""
    return S / np.pi * gamma_L / ((nu - nu0)**2 + gamma_L**2)

def doppler_line(nu, nu0, S, gamma_D):
    """多普勒线型"""
    return S / gamma_D * np.sqrt(np.log(2) / np.pi) * \
           np.exp(-np.log(2) * ((nu - nu0) / gamma_D)**2)

def voigt_line(nu, nu0, S, gamma_L, gamma_D):
    """Voigt线型 (使用wofz函数)"""
    sigma = gamma_D / np.sqrt(2 * np.log(2))
    z = ((nu - nu0) + 1j * gamma_L) / (sigma * np.sqrt(2))
    return S * np.real(wofz(z)) / (sigma * np.sqrt(2 * np.pi))

def case1_line_shapes():
    """案例1：线型对比"""
    print("=" * 60)
    print("案例1：单条谱线线型分析")
    print("=" * 60)
    
    # 谱线参数
    nu0 = 2000.0  # 中心波数 [cm⁻¹]
    S = 1.0  # 谱线强度 [cm⁻² atm⁻¹]
    gamma_L = 0.1  # 洛伦兹半宽度 [cm⁻¹]
    gamma_D = 0.05  # 多普勒半宽度 [cm⁻¹]
    
    # 波数范围
    nu = np.linspace(nu0 - 1, nu0 + 1, 1000)
    
    # 计算各线型
    kappa_L = lorentz_line(nu, nu0, S, gamma_L)
    kappa_D = doppler_line(nu, nu0, S, gamma_D)
    kappa_V = voigt_line(nu, nu0, S, gamma_L, gamma_D)
    
    fig, axes = plt.subplots(2, 2, figsize=(14, 12))
    
    # 图1: 三种线型对比（线性坐标）
    ax1 = axes[0, 0]
    ax1.plot(nu, kappa_L, 'b-', linewidth=2, label='Lorentz')
    ax1.plot(nu, kappa_D, 'r--', linewidth=2, label='Doppler')
    ax1.plot(nu, kappa_V, 'g:', linewidth=2, label='Voigt')
    ax1.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax1.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax1.set_title('线型对比 (线性坐标)', fontsize=14)
    ax1.legend()
    ax1.grid(True, alpha=0.3)
    
    # 图2: 三种线型对比（对数坐标）
    ax2 = axes[0, 1]
    ax2.semilogy(nu, kappa_L, 'b-', linewidth=2, label='Lorentz')
    ax2.semilogy(nu, kappa_D, 'r--', linewidth=2, label='Doppler')
    ax2.semilogy(nu, kappa_V, 'g:', linewidth=2, label='Voigt')
    ax2.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax2.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax2.set_title('线型对比 (对数坐标)', fontsize=14)
    ax2.legend()
    ax2.grid(True, alpha=0.3)
    
    # 图3: 不同温度下的Voigt线型
    ax3 = axes[1, 0]
    T_range = [300, 500, 1000, 1500]
    colors = plt.cm.jet(np.linspace(0, 1, len(T_range)))
    
    for T, color in zip(T_range, colors):
        gamma_D_T = gamma_D * np.sqrt(T / 300)
        gamma_L_T = gamma_L * (300 / T)**0.7
        kappa_V_T = voigt_line(nu, nu0, S, gamma_L_T, gamma_D_T)
        ax3.plot(nu, kappa_V_T, color=color, linewidth=2, label=f'T={T}K')
    
    ax3.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax3.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax3.set_title('不同温度下的Voigt线型', fontsize=14)
    ax3.legend()
    ax3.grid(True, alpha=0.3)
    
    # 图4: 线翼行为对比
    ax4 = axes[1, 1]
    nu_far = np.linspace(nu0 - 10, nu0 + 10, 2000)
    kappa_L_far = lorentz_line(nu_far, nu0, S, gamma_L)
    kappa_D_far = doppler_line(nu_far, nu0, S, gamma_D)
    kappa_V_far = voigt_line(nu_far, nu0, S, gamma_L, gamma_D)
    
    ax4.semilogy(nu_far, kappa_L_far, 'b-', linewidth=2, label='Lorentz')
    ax4.semilogy(nu_far, kappa_D_far, 'r--', linewidth=2, label='Doppler')
    ax4.semilogy(nu_far, kappa_V_far, 'g:', linewidth=2, label='Voigt')
    ax4.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax4.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax4.set_title('线翼行为对比', fontsize=14)
    ax4.legend()
    ax4.grid(True, alpha=0.3)
    ax4.set_ylim(1e-4, 10)
    
    plt.tight_layout()
    plt.savefig('case1_line_shapes.png', dpi=150, bbox_inches='tight')
    plt.close()
    print("案例1完成：已保存 case1_line_shapes.png")

if __name__ == '__main__':
    case1_line_shapes()

案例2：多条谱线叠加

"""
主题027：全谱线-by-线计算仿真程序
包含6个案例和3个GIF动画
"""

import numpy as np
import matplotlib.pyplot as plt
from matplotlib.animation import FuncAnimation, PillowWriter
from scipy.special import wofz
from scipy.integrate import trapezoid
import matplotlib
matplotlib.use('Agg')
import warnings
warnings.filterwarnings('ignore')

# 设置中文字体
plt.rcParams['font.sans-serif'] = ['SimHei', 'DejaVu Sans', 'Arial Unicode MS']
plt.rcParams['axes.unicode_minus'] = False

# 兼容不同numpy版本
try:
    from numpy import trapz
except ImportError:
    from scipy.integrate import trapezoid as trapz

# ==================== 基础函数定义 ====================

def lorentz_line(nu, nu0, S, gamma_L):
    """洛伦兹线型"""
    return S / np.pi * gamma_L / ((nu - nu0)**2 + gamma_L**2)

def doppler_line(nu, nu0, S, gamma_D):
    """多普勒线型"""
    return S / gamma_D * np.sqrt(np.log(2) / np.pi) * \
           np.exp(-np.log(2) * ((nu - nu0) / gamma_D)**2)

def voigt_line(nu, nu0, S, gamma_L, gamma_D):
    """Voigt线型 (使用wofz函数)"""
    sigma = gamma_D / np.sqrt(2 * np.log(2))
    z = ((nu - nu0) + 1j * gamma_L) / (sigma * np.sqrt(2))
    return S * np.real(wofz(z)) / (sigma * np.sqrt(2 * np.pi))

def planck_distribution(nu, T):
    """普朗克分布 (波数形式)"""
    c1 = 1.191042e-5  # mW/(m²·sr·cm⁻⁴)
    c2 = 1.438777  # cm·K
    return c1 * nu**3 / (np.exp(c2 * nu / T) - 1)

def line_strength_temperature(S_ref, T_ref, T, E_low, nu0, Q_ratio):
    """计算谱线强度的温度依赖"""
    c2 = 1.438777  # cm·K
    ratio = Q_ratio * np.exp(-c2 * E_low * (1/T - 1/T_ref)) * \
            (1 - np.exp(-c2 * nu0/T)) / (1 - np.exp(-c2 * nu0/T_ref))
    return S_ref * ratio

# ==================== 案例1：单条谱线线型分析 ====================

def case1_line_shapes():
    """案例1：单条谱线线型分析"""
    print("=" * 60)
    print("案例1：单条谱线线型分析")
    print("=" * 60)
    
    # 谱线参数
    nu0 = 2000.0  # 中心波数 [cm⁻¹]
    S = 1.0  # 谱线强度 [cm⁻² atm⁻¹]
    gamma_L = 0.1  # 洛伦兹半宽度 [cm⁻¹]
    gamma_D = 0.05  # 多普勒半宽度 [cm⁻¹]
    
    # 波数范围
    nu = np.linspace(nu0 - 1, nu0 + 1, 1000)
    
    # 计算各线型
    kappa_L = lorentz_line(nu, nu0, S, gamma_L)
    kappa_D = doppler_line(nu, nu0, S, gamma_D)
    kappa_V = voigt_line(nu, nu0, S, gamma_L, gamma_D)
    
    fig, axes = plt.subplots(2, 2, figsize=(14, 12))
    
    # 图1: 三种线型对比（线性坐标）
    ax1 = axes[0, 0]
    ax1.plot(nu, kappa_L, 'b-', linewidth=2, label='Lorentz')
    ax1.plot(nu, kappa_D, 'r--', linewidth=2, label='Doppler')
    ax1.plot(nu, kappa_V, 'g:', linewidth=2, label='Voigt')
    ax1.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax1.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax1.set_title('线型对比 (线性坐标)', fontsize=14)
    ax1.legend()
    ax1.grid(True, alpha=0.3)
    
    # 图2: 三种线型对比（对数坐标）
    ax2 = axes[0, 1]
    ax2.semilogy(nu, kappa_L, 'b-', linewidth=2, label='Lorentz')
    ax2.semilogy(nu, kappa_D, 'r--', linewidth=2, label='Doppler')
    ax2.semilogy(nu, kappa_V, 'g:', linewidth=2, label='Voigt')
    ax2.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax2.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax2.set_title('线型对比 (对数坐标)', fontsize=14)
    ax2.legend()
    ax2.grid(True, alpha=0.3)
    
    # 图3: 不同温度下的Voigt线型
    ax3 = axes[1, 0]
    T_range = [300, 500, 1000, 1500]
    colors = plt.cm.jet(np.linspace(0, 1, len(T_range)))
    
    for T, color in zip(T_range, colors):
        gamma_D_T = gamma_D * np.sqrt(T / 300)
        gamma_L_T = gamma_L * (300 / T)**0.7
        kappa_V_T = voigt_line(nu, nu0, S, gamma_L_T, gamma_D_T)
        ax3.plot(nu, kappa_V_T, color=color, linewidth=2, label=f'T={T}K')
    
    ax3.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax3.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax3.set_title('不同温度下的Voigt线型', fontsize=14)
    ax3.legend()
    ax3.grid(True, alpha=0.3)
    
    # 图4: 线翼行为对比
    ax4 = axes[1, 1]
    nu_far = np.linspace(nu0 - 10, nu0 + 10, 2000)
    kappa_L_far = lorentz_line(nu_far, nu0, S, gamma_L)
    kappa_D_far = doppler_line(nu_far, nu0, S, gamma_D)
    kappa_V_far = voigt_line(nu_far, nu0, S, gamma_L, gamma_D)
    
    ax4.semilogy(nu_far, kappa_L_far, 'b-', linewidth=2, label='Lorentz')
    ax4.semilogy(nu_far, kappa_D_far, 'r--', linewidth=2, label='Doppler')
    ax4.semilogy(nu_far, kappa_V_far, 'g:', linewidth=2, label='Voigt')
    ax4.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax4.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax4.set_title('线翼行为对比', fontsize=14)
    ax4.legend()
    ax4.grid(True, alpha=0.3)
    ax4.set_ylim(1e-4, 10)
    
    plt.tight_layout()
    plt.savefig('case1_line_shapes.png', dpi=150, bbox_inches='tight')
    plt.close()
    print("案例1完成：已保存 case1_line_shapes.png")

# ==================== 案例2：多条谱线叠加 ====================

def case2_multiple_lines():
    """案例2：多条谱线叠加模拟"""
    print("=" * 60)
    print("案例2：多条谱线叠加模拟")
    print("=" * 60)
    
    # 模拟CO2 4.3μm带的多条谱线 (简化)
    # 生成随机谱线参数
    np.random.seed(42)
    n_lines = 50
    nu_center = 2300  # 中心波数
    nu_range = 100  # 范围
    
    line_data = []
    for i in range(n_lines):
        nu0 = nu_center + np.random.uniform(-nu_range, nu_range)
        S = np.random.uniform(0.1, 2.0) * np.exp(-0.5 * ((nu0 - nu_center) / 50)**2)
        gamma_L = 0.08
        gamma_D = 0.03
        line_data.append({'nu0': nu0, 'S': S, 'gamma_L': gamma_L, 'gamma_D': gamma_D})
    
    # 波数网格
    nu = np.linspace(nu_center - nu_range, nu_center + nu_range, 2000)
    
    # 计算总吸收系数
    kappa_total = np.zeros_like(nu)
    for line in line_data:
        kappa_total += voigt_line(nu, line['nu0'], line['S'], 
                                  line['gamma_L'], line['gamma_D'])
    
    fig, axes = plt.subplots(2, 2, figsize=(14, 12))
    
    # 图1: 总吸收系数
    ax1 = axes[0, 0]
    ax1.fill_between(nu, kappa_total, alpha=0.5, color='blue')
    ax1.plot(nu, kappa_total, 'b-', linewidth=1)
    ax1.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax1.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax1.set_title('多条谱线叠加的吸收系数', fontsize=14)
    ax1.grid(True, alpha=0.3)
    
    # 图2: 对数坐标
    ax2 = axes[0, 1]
    ax2.semilogy(nu, kappa_total, 'b-', linewidth=1)
    ax2.fill_between(nu, kappa_total, alpha=0.5, color='blue')
    ax2.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax2.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax2.set_title('吸收系数 (对数坐标)', fontsize=14)
    ax2.grid(True, alpha=0.3)
    
    # 图3: 逐条谱线贡献
    ax3 = axes[1, 0]
    colors = plt.cm.rainbow(np.linspace(0, 1, min(10, n_lines)))
    for i, (line, color) in enumerate(zip(line_data[:10], colors)):
        kappa_single = voigt_line(nu, line['nu0'], line['S'],
                                  line['gamma_L'], line['gamma_D'])
        ax3.plot(nu, kappa_single, color=color, linewidth=1, alpha=0.7, 
                label=f'Line {i+1}')
    ax3.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax3.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax3.set_title('前10条谱线的单独贡献', fontsize=14)
    ax3.legend(fontsize=8)
    ax3.grid(True, alpha=0.3)
    
    # 图4: 谱线强度分布
    ax4 = axes[1, 1]
    nu0_list = [line['nu0'] for line in line_data]
    S_list = [line['S'] for line in line_data]
    ax4.scatter(nu0_list, S_list, c='blue', alpha=0.6, s=50)
    ax4.set_xlabel('谱线中心波数 ν₀ [cm⁻¹]', fontsize=12)
    ax4.set_ylabel('谱线强度 S [cm⁻² atm⁻¹]', fontsize=12)
    ax4.set_title('谱线强度分布', fontsize=14)
    ax4.grid(True, alpha=0.3)
    
    plt.tight_layout()
    plt.savefig('case2_multiple_lines.png', dpi=150, bbox_inches='tight')
    plt.close()
    print("案例2完成：已保存 case2_multiple_lines.png")

# ==================== 案例3：光谱分辨率的影响 ====================

def case3_resolution_effect():
    """案例3：光谱分辨率对计算结果的影响"""
    print("=" * 60)
    print("案例3：光谱分辨率的影响")
    print("=" * 60)
    
    # 单条谱线
    nu0 = 2000.0
    S = 1.0
    gamma_L = 0.1
    gamma_D = 0.05
    
    # 不同分辨率
    resolutions = [0.001, 0.01, 0.05, 0.1]  # cm⁻¹
    
    fig, axes = plt.subplots(2, 2, figsize=(14, 12))
    
    # 图1: 不同分辨率的线型
    ax1 = axes[0, 0]
    colors = plt.cm.jet(np.linspace(0, 1, len(resolutions)))
    
    for res, color in zip(resolutions, colors):
        nu = np.arange(nu0 - 1, nu0 + 1, res)
        kappa = voigt_line(nu, nu0, S, gamma_L, gamma_D)
        ax1.plot(nu, kappa, color=color, linewidth=2, label=f'Δν={res} cm⁻¹')
    
    # 高分辨率参考
    nu_fine = np.linspace(nu0 - 1, nu0 + 1, 10000)
    kappa_fine = voigt_line(nu_fine, nu0, S, gamma_L, gamma_D)
    ax1.plot(nu_fine, kappa_fine, 'k--', linewidth=1, alpha=0.5, label='参考')
    
    ax1.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax1.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax1.set_title('不同分辨率的线型', fontsize=14)
    ax1.legend()
    ax1.grid(True, alpha=0.3)
    
    # 图2: 积分误差
    ax2 = axes[0, 1]
    # 计算等效线宽
    L = 1.0  # 行程长度
    W_exact = trapz(1 - np.exp(-kappa_fine * L), nu_fine)
    
    errors = []
    for res in resolutions:
        nu = np.arange(nu0 - 1, nu0 + 1, res)
        kappa = voigt_line(nu, nu0, S, gamma_L, gamma_D)
        W = trapz(1 - np.exp(-kappa * L), nu)
        error = abs(W - W_exact) / W_exact * 100
        errors.append(error)
    
    ax2.semilogy(resolutions, errors, 'bo-', linewidth=2, markersize=8)
    ax2.set_xlabel('波数分辨率 Δν [cm⁻¹]', fontsize=12)
    ax2.set_ylabel('等效线宽误差 [%]', fontsize=12)
    ax2.set_title('分辨率对积分精度的影响', fontsize=14)
    ax2.grid(True, alpha=0.3)
    
    # 图3: 计算时间对比 (模拟)
    ax3 = axes[1, 0]
    n_points = [2000 / res for res in resolutions]
    times = [n / 1000 for n in n_points]  # 模拟时间
    
    ax3.bar(range(len(resolutions)), times, color=colors, alpha=0.7, edgecolor='black')
    ax3.set_xticks(range(len(resolutions)))
    ax3.set_xticklabels([f'{res}' for res in resolutions])
    ax3.set_xlabel('波数分辨率 Δν [cm⁻¹]', fontsize=12)
    ax3.set_ylabel('相对计算时间', fontsize=12)
    ax3.set_title('分辨率对计算时间的影响', fontsize=14)
    ax3.grid(True, alpha=0.3, axis='y')
    
    # 图4: 精度-效率权衡
    ax4 = axes[1, 1]
    ax4.scatter(times, errors, s=200, c=colors, alpha=0.7, edgecolors='black')
    for i, res in enumerate(resolutions):
        ax4.annotate(f'{res}', (times[i], errors[i]), 
                    xytext=(10, 10), textcoords='offset points', fontsize=10)
    ax4.set_xlabel('相对计算时间', fontsize=12)
    ax4.set_ylabel('误差 [%]', fontsize=12)
    ax4.set_title('精度-效率权衡', fontsize=14)
    ax4.set_yscale('log')
    ax4.grid(True, alpha=0.3)
    
    plt.tight_layout()
    plt.savefig('case3_resolution_effect.png', dpi=150, bbox_inches='tight')
    plt.close()
    print("案例3完成：已保存 case3_resolution_effect.png")

# ==================== 案例4：温度对光谱的影响 ====================

def case4_temperature_effect():
    """案例4：温度对光谱的影响"""
    print("=" * 60)
    print("案例4：温度对光谱的影响")
    print("=" * 60)
    
    # 模拟多条谱线
    np.random.seed(42)
    n_lines = 30
    nu_center = 2000
    
    # 基准谱线参数 (300K)
    line_data = []
    for i in range(n_lines):
        nu0 = nu0 = nu_center + np.random.uniform(-100, 100)
        E_low = np.random.uniform(0, 2000)  # 低能级能量
        S_300 = np.random.uniform(0.1, 1.0)
        line_data.append({'nu0': nu0, 'E_low': E_low, 'S_300': S_300})
    
    T_range = [300, 600, 1000, 1500]
    nu = np.linspace(nu_center - 150, nu_center + 150, 3000)
    
    fig, axes = plt.subplots(2, 2, figsize=(14, 12))
    
    # 图1: 不同温度下的吸收系数
    ax1 = axes[0, 0]
    colors = plt.cm.jet(np.linspace(0, 1, len(T_range)))
    
    for T, color in zip(T_range, colors):
        kappa = np.zeros_like(nu)
        for line in line_data:
            S_T = line_strength_temperature(line['S_300'], 300, T, 
                                           line['E_low'], line['nu0'], 1.0)
            gamma_L = 0.1 * (300 / T)**0.7
            gamma_D = 0.03 * np.sqrt(T / 300)
            kappa += voigt_line(nu, line['nu0'], S_T, gamma_L, gamma_D)
        
        ax1.plot(nu, kappa, color=color, linewidth=1.5, label=f'T={T}K')
    
    ax1.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax1.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax1.set_title('不同温度下的吸收系数', fontsize=14)
    ax1.legend()
    ax1.grid(True, alpha=0.3)
    
    # 图2: 总吸收对比
    ax2 = axes[0, 1]
    total_absorption = []
    for T in T_range:
        kappa = np.zeros_like(nu)
        for line in line_data:
            S_T = line_strength_temperature(line['S_300'], 300, T, 
                                           line['E_low'], line['nu0'], 1.0)
            gamma_L = 0.1 * (300 / T)**0.7
            gamma_D = 0.03 * np.sqrt(T / 300)
            kappa += voigt_line(nu, line['nu0'], S_T, gamma_L, gamma_D)
        
        # 计算总吸收 (简化)
        A = trapz(kappa, nu)
        total_absorption.append(A)
    
    ax2.plot(T_range, total_absorption, 'ro-', linewidth=2, markersize=8)
    ax2.set_xlabel('温度 T [K]', fontsize=12)
    ax2.set_ylabel('总吸收 ∫κ dν [cm⁻¹ atm⁻¹]', fontsize=12)
    ax2.set_title('总吸收随温度的变化', fontsize=14)
    ax2.grid(True, alpha=0.3)
    
    # 图3: 谱线强度分布变化
    ax3 = axes[1, 0]
    T_plot = [300, 1000, 1500]
    colors_T = ['blue', 'green', 'red']
    
    for T, color in zip(T_plot, colors_T):
        S_values = []
        for line in line_data:
            S_T = line_strength_temperature(line['S_300'], 300, T, 
                                           line['E_low'], line['nu0'], 1.0)
            S_values.append(S_T)
        
        nu0_list = [line['nu0'] for line in line_data]
        ax3.scatter(nu0_list, S_values, c=color, alpha=0.6, s=50, label=f'T={T}K')
    
    ax3.set_xlabel('谱线中心波数 ν₀ [cm⁻¹]', fontsize=12)
    ax3.set_ylabel('谱线强度 S [cm⁻² atm⁻¹]', fontsize=12)
    ax3.set_title('谱线强度随温度的变化', fontsize=14)
    ax3.legend()
    ax3.grid(True, alpha=0.3)
    
    # 图4: 发射率计算
    ax4 = axes[1, 1]
    L = 1.0  # 行程长度
    epsilon = []
    
    for T in T_range:
        kappa = np.zeros_like(nu)
        for line in line_data:
            S_T = line_strength_temperature(line['S_300'], 300, T, 
                                           line['E_low'], line['nu0'], 1.0)
            gamma_L = 0.1 * (300 / T)**0.7
            gamma_D = 0.03 * np.sqrt(T / 300)
            kappa += voigt_line(nu, line['nu0'], S_T, gamma_L, gamma_D)
        
        # 计算发射率 (简化)
        tau = np.exp(-kappa * L)
        eps = 1 - np.mean(tau)
        epsilon.append(eps)
    
    ax4.plot(T_range, epsilon, 'bs-', linewidth=2, markersize=8)
    ax4.set_xlabel('温度 T [K]', fontsize=12)
    ax4.set_ylabel('发射率 ε', fontsize=12)
    ax4.set_title('发射率随温度的变化', fontsize=14)
    ax4.grid(True, alpha=0.3)
    ax4.set_ylim(0, 1)
    
    plt.tight_layout()
    plt.savefig('case4_temperature_effect.png', dpi=150, bbox_inches='tight')
    plt.close()
    print("案例4完成：已保存 case4_temperature_effect.png")

# ==================== 案例5：LBL vs 简化模型 ====================

def case5_lbl_vs_models():
    """案例5：LBL计算与简化模型对比"""
    print("=" * 60)
    print("案例5：LBL与简化模型对比")
    print("=" * 60)
    
    # 生成模拟谱线 (作为"LBL"基准)
    np.random.seed(42)
    n_lines = 100
    nu_center = 2000
    nu_range = 200
    
    line_data = []
    for i in range(n_lines):
        nu0 = nu_center + np.random.uniform(-nu_range, nu_range)
        S = np.random.uniform(0.05, 0.5) * np.exp(-0.5 * ((nu0 - nu_center) / 100)**2)
        gamma_L = 0.08
        gamma_D = 0.03
        line_data.append({'nu0': nu0, 'S': S, 'gamma_L': gamma_L, 'gamma_D': gamma_D})
    
    # 高分辨率计算 (LBL)
    nu_fine = np.linspace(nu_center - nu_range, nu_center + nu_range, 5000)
    kappa_lbl = np.zeros_like(nu_fine)
    for line in line_data:
        kappa_lbl += voigt_line(nu_fine, line['nu0'], line['S'],
                                line['gamma_L'], line['gamma_D'])
    
    # 低分辨率近似 (模拟窄谱带模型)
    nu_coarse = np.linspace(nu_center - nu_range, nu_center + nu_range, 50)
    kappa_coarse = np.zeros_like(nu_coarse)
    for i, nu_c in enumerate(nu_coarse):
        mask = (nu_fine >= nu_c - 2) & (nu_fine < nu_c + 2)
        if np.any(mask):
            kappa_coarse[i] = np.mean(kappa_lbl[mask])
    
    # 宽谱带近似 (模拟WBM)
    kappa_wide = np.ones_like(nu_fine) * np.mean(kappa_lbl)
    
    fig, axes = plt.subplots(2, 2, figsize=(14, 12))
    
    # 图1: 吸收系数对比
    ax1 = axes[0, 0]
    ax1.plot(nu_fine, kappa_lbl, 'b-', linewidth=1, label='LBL (高分辨率)', alpha=0.8)
    ax1.plot(nu_coarse, kappa_coarse, 'ro-', linewidth=2, markersize=4, 
             label='窄谱带近似', alpha=0.7)
    ax1.plot(nu_fine, kappa_wide, 'g--', linewidth=2, label='宽谱带近似', alpha=0.7)
    ax1.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax1.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax1.set_title('不同模型的吸收系数对比', fontsize=14)
    ax1.legend()
    ax1.grid(True, alpha=0.3)
    
    # 图2: 透射率对比
    ax2 = axes[0, 1]
    L = 1.0
    tau_lbl = np.exp(-kappa_lbl * L)
    tau_coarse = np.exp(-kappa_coarse * L)
    tau_wide = np.exp(-kappa_wide * L)
    
    ax2.plot(nu_fine, tau_lbl, 'b-', linewidth=1, label='LBL', alpha=0.8)
    ax2.plot(nu_coarse, tau_coarse, 'ro-', linewidth=2, markersize=4, 
             label='窄谱带', alpha=0.7)
    ax2.plot(nu_fine, tau_wide, 'g--', linewidth=2, label='宽谱带', alpha=0.7)
    ax2.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax2.set_ylabel('透射率 τ', fontsize=12)
    ax2.set_title('透射率对比', fontsize=14)
    ax2.legend()
    ax2.grid(True, alpha=0.3)
    
    # 图3: 发射率误差分析
    ax3 = axes[1, 0]
    L_range = np.logspace(-2, 1, 50)
    
    eps_lbl = []
    eps_coarse = []
    eps_wide = []
    
    for L_val in L_range:
        tau_l = np.exp(-kappa_lbl * L_val)
        tau_c = np.exp(-kappa_coarse * L_val)
        tau_w = np.exp(-kappa_wide * L_val)
        
        # 简化发射率计算
        eps_lbl.append(1 - np.mean(tau_l))
        eps_coarse.append(1 - np.mean(tau_c))
        eps_wide.append(1 - np.mean(tau_w))
    
    error_coarse = np.abs(np.array(eps_coarse) - np.array(eps_lbl)) / np.array(eps_lbl) * 100
    error_wide = np.abs(np.array(eps_wide) - np.array(eps_lbl)) / np.array(eps_lbl) * 100
    
    ax3.semilogx(L_range, error_coarse, 'r-', linewidth=2, label='窄谱带误差')
    ax3.semilogx(L_range, error_wide, 'g--', linewidth=2, label='宽谱带误差')
    ax3.set_xlabel('行程长度 L [m]', fontsize=12)
    ax3.set_ylabel('发射率相对误差 [%]', fontsize=12)
    ax3.set_title('简化模型误差分析', fontsize=14)
    ax3.legend()
    ax3.grid(True, alpha=0.3)
    
    # 图4: 计算成本对比
    ax4 = axes[1, 1]
    models = ['LBL', '窄谱带', '宽谱带', 'WSGGM']
    costs = [1000, 50, 5, 1]
    errors = [0, np.mean(error_coarse), np.mean(error_wide), 30]
    colors = ['blue', 'orange', 'green', 'red']
    
    scatter = ax4.scatter(costs, errors, s=300, c=colors, alpha=0.7, edgecolors='black')
    for i, model in enumerate(models):
        ax4.annotate(model, (costs[i], errors[i]), 
                    xytext=(10, 10), textcoords='offset points', fontsize=11)
    
    ax4.set_xlabel('相对计算成本', fontsize=12)
    ax4.set_ylabel('平均误差 [%]', fontsize=12)
    ax4.set_title('精度-成本权衡', fontsize=14)
    ax4.set_xscale('log')
    ax4.grid(True, alpha=0.3)
    
    plt.tight_layout()
    plt.savefig('case5_lbl_vs_models.png', dpi=150, bbox_inches='tight')
    plt.close()
    print("案例5完成：已保存 case5_lbl_vs_models.png")

# ==================== 案例6：高分辨率光谱分析 ====================

def case6_high_resolution():
    """案例6：高分辨率光谱分析"""
    print("=" * 60)
    print("案例6：高分辨率光谱分析")
    print("=" * 60)
    
    # 模拟高分辨率光谱 (类似HITRAN数据)
    np.random.seed(123)
    nu_center = 2000
    
    # 生成密集的谱线
    n_lines = 500
    line_data = []
    for i in range(n_lines):
        nu0 = nu_center + np.random.normal(0, 150)
        S = np.random.lognormal(-3, 1.5)
        E_low = np.random.uniform(0, 3000)
        gamma_L = 0.08
        gamma_D = 0.03
        line_data.append({'nu0': nu0, 'S': S, 'E_low': E_low, 
                         'gamma_L': gamma_L, 'gamma_D': gamma_D})
    
    # 不同分辨率
    resolutions = {
        '高分辨率 (0.001 cm⁻¹)': 0.001,
        '中分辨率 (0.01 cm⁻¹)': 0.01,
        '低分辨率 (0.1 cm⁻¹)': 0.1
    }
    
    T = 1200  # K
    
    fig, axes = plt.subplots(2, 2, figsize=(14, 12))
    
    # 图1: 不同分辨率的光谱
    ax1 = axes[0, 0]
    colors = ['blue', 'green', 'red']
    
    for (name, res), color in zip(resolutions.items(), colors):
        nu = np.arange(nu_center - 100, nu_center + 100, res)
        kappa = np.zeros_like(nu)
        
        for line in line_data:
            S_T = line_strength_temperature(line['S'], 300, T, 
                                           line['E_low'], line['nu0'], 1.0)
            gamma_L_T = line['gamma_L'] * (300 / T)**0.7
            gamma_D_T = line['gamma_D'] * np.sqrt(T / 300)
            kappa += voigt_line(nu, line['nu0'], S_T, gamma_L_T, gamma_D_T)
        
        ax1.plot(nu, kappa, color=color, linewidth=0.8 if res < 0.01 else 1.5, 
                label=name, alpha=0.8)
    
    ax1.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax1.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
    ax1.set_title('不同分辨率的光谱', fontsize=14)
    ax1.legend()
    ax1.grid(True, alpha=0.3)
    ax1.set_xlim(nu_center - 50, nu_center + 50)
    
    # 图2: 谱线密度分布
    ax2 = axes[0, 1]
    nu0_list = [line['nu0'] for line in line_data]
    ax2.hist(nu0_list, bins=50, color='blue', alpha=0.7, edgecolor='black')
    ax2.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
    ax2.set_ylabel('谱线数量', fontsize=12)
    ax2.set_title('谱线密度分布', fontsize=14)
    ax2.grid(True, alpha=0.3)
    
    # 图3: 谱线强度分布
    ax3 = axes[1, 0]
    S_list = [line['S'] for line in line_data]
    ax3.hist(np.log10(S_list), bins=30, color='green', alpha=0.7, edgecolor='black')
    ax3.set_xlabel('log₁₀(谱线强度)', fontsize=12)
    ax3.set_ylabel('谱线数量', fontsize=12)
    ax3.set_title('谱线强度分布', fontsize=14)
    ax3.grid(True, alpha=0.3)
    
    # 图4: 累积吸收贡献
    ax4 = axes[1, 1]
    
    # 计算高分辨率光谱
    nu_high = np.linspace(nu_center - 200, nu_center + 200, 10000)
    kappa_high = np.zeros_like(nu_high)
    
    for line in line_data:
        S_T = line_strength_temperature(line['S'], 300, T, 
                                       line['E_low'], line['nu0'], 1.0)
        gamma_L_T = line['gamma_L'] * (300 / T)**0.7
        gamma_D_T = line['gamma_D'] * np.sqrt(T / 300)
        kappa_high += voigt_line(nu_high, line['nu0'], S_T, gamma_L_T, gamma_D_T)
    
    # 计算累积分布
    kappa_sorted = np.sort(kappa_high)[::-1]
    cumulative = np.cumsum(kappa_sorted) / np.sum(kappa_sorted) * 100
    x_percent = np.linspace(0, 100, len(cumulative))
    
    ax4.plot(x_percent, cumulative, 'b-', linewidth=2)
    ax4.axhline(y=90, color='r', linestyle='--', alpha=0.7, label='90%吸收')
    ax4.axvline(x=10, color='g', linestyle='--', alpha=0.7, label='10%波数范围')
    ax4.set_xlabel('波数范围 [%]', fontsize=12)
    ax4.set_ylabel('累积吸收 [%]', fontsize=12)
    ax4.set_title('累积吸收贡献', fontsize=14)
    ax4.legend()
    ax4.grid(True, alpha=0.3)
    
    plt.tight_layout()
    plt.savefig('case6_high_resolution.png', dpi=150, bbox_inches='tight')
    plt.close()
    print("案例6完成：已保存 case6_high_resolution.png")

# ==================== GIF动画生成 ====================

def create_gif_line_shape_evolution():
    """创建GIF动画1：线型随温度演化"""
    print("正在生成GIF动画1：线型演化...")
    
    nu0 = 2000.0
    S = 1.0
    gamma_L_ref = 0.1
    gamma_D_ref = 0.05
    
    T_range = np.linspace(300, 1500, 40)
    nu = np.linspace(nu0 - 1, nu0 + 1, 500)
    
    fig, ax = plt.subplots(figsize=(10, 6))
    
    def update(frame):
        ax.clear()
        
        T = T_range[frame]
        gamma_L = gamma_L_ref * (300 / T)**0.7
        gamma_D = gamma_D_ref * np.sqrt(T / 300)
        
        kappa_L = lorentz_line(nu, nu0, S, gamma_L)
        kappa_D = doppler_line(nu, nu0, S, gamma_D)
        kappa_V = voigt_line(nu, nu0, S, gamma_L, gamma_D)
        
        ax.plot(nu, kappa_L, 'b-', linewidth=2, label='Lorentz', alpha=0.7)
        ax.plot(nu, kappa_D, 'r--', linewidth=2, label='Doppler', alpha=0.7)
        ax.plot(nu, kappa_V, 'g-', linewidth=2.5, label='Voigt')
        
        ax.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
        ax.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
        ax.set_title(f'线型演化 (T = {T:.0f} K)', fontsize=14)
        ax.legend()
        ax.grid(True, alpha=0.3)
        ax.set_ylim(0, 4)
        
        return ax,
    
    ani = FuncAnimation(fig, update, frames=len(T_range), interval=150, blit=False)
    writer = PillowWriter(fps=6)
    ani.save('gif1_line_shape_evolution.gif', writer=writer)
    plt.close()
    print("GIF动画1完成：已保存 gif1_line_shape_evolution.gif")

def create_gif_spectrum_resolution():
    """创建GIF动画2：光谱分辨率变化"""
    print("正在生成GIF动画2：光谱分辨率变化...")
    
    # 生成几条谱线
    np.random.seed(42)
    line_data = [
        {'nu0': 1998, 'S': 1.0, 'gamma_L': 0.1, 'gamma_D': 0.05},
        {'nu0': 2000, 'S': 1.5, 'gamma_L': 0.1, 'gamma_D': 0.05},
        {'nu0': 2002, 'S': 0.8, 'gamma_L': 0.1, 'gamma_D': 0.05},
    ]
    
    resolutions = np.logspace(-3, -1, 30)  # 0.001 to 0.1
    
    fig, ax = plt.subplots(figsize=(10, 6))
    
    def update(frame):
        ax.clear()
        
        res = resolutions[frame]
        nu = np.arange(1995, 2005, res)
        
        kappa = np.zeros_like(nu)
        for line in line_data:
            kappa += voigt_line(nu, line['nu0'], line['S'], 
                               line['gamma_L'], line['gamma_D'])
        
        ax.fill_between(nu, kappa, alpha=0.5, color='blue')
        ax.plot(nu, kappa, 'b-', linewidth=1.5)
        ax.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
        ax.set_ylabel('吸收系数 κ [cm⁻² atm⁻¹]', fontsize=12)
        ax.set_title(f'分辨率变化 (Δν = {res:.4f} cm⁻¹)', fontsize=14)
        ax.grid(True, alpha=0.3)
        ax.set_xlim(1995, 2005)
        
        return ax,
    
    ani = FuncAnimation(fig, update, frames=len(resolutions), interval=200, blit=False)
    writer = PillowWriter(fps=5)
    ani.save('gif2_spectrum_resolution.gif', writer=writer)
    plt.close()
    print("GIF动画2完成：已保存 gif2_spectrum_resolution.gif")

def create_gif_transmittance_evolution():
    """创建GIF动画3：透射率随路径演化"""
    print("正在生成GIF动画3：透射率演化...")
    
    # 生成模拟光谱
    np.random.seed(42)
    n_lines = 50
    nu_center = 2000
    line_data = []
    for i in range(n_lines):
        nu0 = nu_center + np.random.uniform(-50, 50)
        S = np.random.uniform(0.1, 0.5)
        line_data.append({'nu0': nu0, 'S': S, 'gamma_L': 0.08, 'gamma_D': 0.03})
    
    nu = np.linspace(nu_center - 60, nu_center + 60, 1000)
    kappa = np.zeros_like(nu)
    for line in line_data:
        kappa += voigt_line(nu, line['nu0'], line['S'], 
                           line['gamma_L'], line['gamma_D'])
    
    L_range = np.logspace(-2, 1, 40)
    
    fig, ax = plt.subplots(figsize=(10, 6))
    
    def update(frame):
        ax.clear()
        
        L = L_range[frame]
        tau = np.exp(-kappa * L)
        
        ax.fill_between(nu, tau, alpha=0.5, color='green')
        ax.plot(nu, tau, 'g-', linewidth=1.5)
        ax.set_xlabel('波数 ν [cm⁻¹]', fontsize=12)
        ax.set_ylabel('透射率 τ', fontsize=12)
        ax.set_title(f'透射率演化 (L = {L:.3f} m)', fontsize=14)
        ax.grid(True, alpha=0.3)
        ax.set_ylim(0, 1)
        
        return ax,
    
    ani = FuncAnimation(fig, update, frames=len(L_range), interval=150, blit=False)
    writer = PillowWriter(fps=6)
    ani.save('gif3_transmittance_evolution.gif', writer=writer)
    plt.close()
    print("GIF动画3完成：已保存 gif3_transmittance_evolution.gif")

# ==================== 主程序 ====================

if __name__ == '__main__':
    print("=" * 70)
    print("主题027：全谱线-by-线计算仿真程序")
    print("=" * 70)
    
    # 运行所有案例
    case1_line_shapes()
    case2_multiple_lines()
    case3_resolution_effect()
    case4_temperature_effect()
    case5_lbl_vs_models()
    case6_high_resolution()
    
    print("\n" + "=" * 70)
    print("所有静态图片生成完成！")
    print("=" * 70)
    
    # 生成GIF动画
    print("\n开始生成GIF动画...")
    create_gif_line_shape_evolution()
    create_gif_spectrum_resolution()
    create_gif_transmittance_evolution()
    
    print("\n" + "=" * 70)
    print("所有仿真案例和GIF动画生成完成！")
    print("=" * 70)
    
    print("\n生成的文件列表：")
    print("  - case1_line_shapes.png")
    print("  - case2_multiple_lines.png")
    print("  - case3_resolution_effect.png")
    print("  - case4_temperature_effect.png")
    print("  - case5_lbl_vs_models.png")
    print("  - case6_high_resolution.png")
    print("  - gif1_line_shape_evolution.gif")
    print("  - gif2_spectrum_resolution.gif")
    print("  - gif3_transmittance_evolution.gif")
    print("=" * 70)