AstroResearch/skills/plotting/SKILL.md

---
name: plotting
description: 科研绘图规范 —— 天文学常用图表的 Python matplotlib 绘制模板
context: fork
allowed-tools:
  - bash
  - save_note
  - write
  - read
---

# 科研绘图规范

## 环境要求

- Python 3.8+
- matplotlib >= 3.7
- numpy
- 可选：seaborn（统计图）、astropy（天文学单位和坐标）

安装：
```bash
pip install matplotlib numpy seaborn astropy --quiet
```

## 输出规范

| 目标 | 分辨率 | 格式 | 说明 |
|------|--------|------|------|
| 论文投稿 | 300+ dpi | PDF/SVG（矢量优先） | APS/A&A/MNRAS 标准 |
| 演示文稿 | 150 dpi | PNG | 宽度 ≥ 1200px |
| 快速预览 | 100 dpi | PNG | 屏幕查看 |

## 通用样式模板

```python
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import numpy as np

# 天文学论文标准样式
plt.rcParams.update({
    'font.family': 'serif',
    'font.size': 12,
    'axes.labelsize': 14,
    'axes.titlesize': 14,
    'xtick.labelsize': 11,
    'ytick.labelsize': 11,
    'legend.fontsize': 11,
    'figure.figsize': (8, 6),
    'figure.dpi': 300,
    'savefig.dpi': 300,
    'savefig.bbox': 'tight',
    'savefig.pad_inches': 0.05,
    'lines.linewidth': 1.5,
    'axes.linewidth': 1.0,
    'xtick.major.width': 0.8,
    'ytick.major.width': 0.8,
    'xtick.minor.width': 0.5,
    'ytick.minor.width': 0.5,
    'xtick.direction': 'in',
    'ytick.direction': 'in',
    'xtick.top': True,
    'ytick.right': True,
})
```

## 图表模板

### 1. 光谱图

```python
def plot_spectrum(wavelength, flux, title="Spectrum", xlabel=r"Wavelength ($\AA$)", ylabel=r"Flux (erg/s/cm$^2$/$\AA$)", save_path=None):
    fig, ax = plt.subplots()
    ax.plot(wavelength, flux, 'k-', linewidth=0.8)
    ax.set_xlabel(xlabel)
    ax.set_ylabel(ylabel)
    ax.set_title(title)
    ax.minorticks_on()
    if save_path:
        fig.savefig(save_path)
        plt.close(fig)
    else:
        plt.show()
    return fig, ax
```

### 2. 光变曲线

```python
def plot_light_curve(time, flux, time_err=None, flux_err=None, title="Light Curve", xlabel="Time (MJD)", ylabel="Flux", save_path=None):
    fig, ax = plt.subplots()
    if flux_err is not None:
        ax.errorbar(time, flux, yerr=flux_err, xerr=time_err, fmt='o', markersize=3, capsize=2, color='black', ecolor='gray')
    else:
        ax.plot(time, flux, 'ko', markersize=3)
    ax.set_xlabel(xlabel)
    ax.set_ylabel(ylabel)
    ax.set_title(title)
    ax.minorticks_on()
    if save_path:
        fig.savefig(save_path)
        plt.close(fig)
    else:
        plt.show()
    return fig, ax
```

### 3. 折叠光变曲线

```python
def plot_folded_lc(phase, flux, flux_err=None, period=None, title="Folded Light Curve", save_path=None):
    fig, ax = plt.subplots()
    if flux_err is not None:
        ax.errorbar(phase, flux, yerr=flux_err, fmt='o', markersize=2, capsize=1, color='black', ecolor='gray', alpha=0.7)
    else:
        ax.plot(phase, flux, 'ko', markersize=2, alpha=0.7)
    ax.set_xlabel("Phase")
    ax.set_ylabel("Flux")
    if period:
        ax.set_title(f"{title} (P = {period:.4f} d)")
    else:
        ax.set_title(title)
    ax.set_xlim(0, 1)
    ax.minorticks_on()
    if save_path:
        fig.savefig(save_path)
        plt.close(fig)
    else:
        plt.show()
    return fig, ax
```

### 4. 赫罗图（CMD）

```python
def plot_hr_diagram(bp_rp, abs_g, title="HR Diagram", save_path=None, color_by_density=False):
    fig, ax = plt.subplots()
    if color_by_density and len(bp_rp) > 100:
        from scipy.stats import gaussian_kde
        xy = np.vstack([bp_rp, abs_g])
        z = gaussian_kde(xy)(xy)
        idx = z.argsort()
        ax.scatter(bp_rp[idx], abs_g[idx], c=z[idx], s=5, cmap='viridis_r', edgecolors='none')
        cb = fig.colorbar(ax.collections[0], ax=ax, pad=0.02)
        cb.set_label("Stellar density")
    else:
        ax.plot(bp_rp, abs_g, 'k.', markersize=1, alpha=0.5)
    ax.set_xlabel(r"$G_{BP} - G_{RP}$ (mag)")
    ax.set_ylabel(r"$M_G$ (mag)")
    ax.set_title(title)
    ax.invert_yaxis()
    ax.minorticks_on()
    if save_path:
        fig.savefig(save_path)
        plt.close(fig)
    else:
        plt.show()
    return fig, ax
```

### 5. SED 图

```python
def plot_sed(wavelength_angstrom, flux_obs, flux_err=None, flux_model=None, title="SED", save_path=None):
    fig, ax = plt.subplots()
    ax.scatter(wavelength_angstrom, flux_obs, c='black', s=30, zorder=5, label='Observed')
    if flux_err is not None:
        ax.errorbar(wavelength_angstrom, flux_obs, yerr=flux_err, fmt='none', ecolor='gray', capsize=3)
    if flux_model is not None:
        model_wave, model_flux = zip(*flux_model) if isinstance(flux_model, list) else (flux_model[0], flux_model[1])
        ax.plot(model_wave, model_flux, 'r-', linewidth=1.5, label='Model', alpha=0.8)
    ax.set_xscale('log')
    ax.set_yscale('log')
    ax.set_xlabel(r"Wavelength ($\AA$)")
    ax.set_ylabel(r"Flux")
    ax.set_title(title)
    ax.legend()
    ax.minorticks_on()
    if save_path:
        fig.savefig(save_path)
        plt.close(fig)
    else:
        plt.show()
    return fig, ax
```

### 6. [α/Fe] vs [Fe/H] 图

```python
def plot_abundance(feh, alpha_feh, labels=None, title="[α/Fe] vs [Fe/H]", save_path=None):
    fig, ax = plt.subplots()
    if labels is not None:
        from matplotlib.colors import ListedColormap
        colors = ['#3498db', '#e74c3c', '#2ecc71']
        cmap = ListedColormap(colors[:len(set(labels))])
        unique = sorted(set(labels))
        for i, lab in enumerate(unique):
            mask = np.array(labels) == lab
            ax.scatter(np.array(feh)[mask], np.array(alpha_feh)[mask], s=10, alpha=0.7, label=lab, color=colors[i % len(colors)])
        ax.legend()
    else:
        ax.scatter(feh, alpha_feh, s=10, alpha=0.7, c='black')
    ax.set_xlabel("[Fe/H] (dex)")
    ax.set_ylabel(r"[$\alpha$/Fe] (dex)")
    ax.set_title(title)
    ax.axhline(y=0.25, color='gray', linestyle='--', linewidth=0.8, alpha=0.5)
    ax.minorticks_on()
    if save_path:
        fig.savefig(save_path)
        plt.close(fig)
    else:
        plt.show()
    return fig, ax
```

### 7. 周期图

```python
def plot_periodogram(periods, powers, best_period=None, title="Periodogram", save_path=None):
    fig, ax = plt.subplots()
    ax.plot(periods, powers, 'k-', linewidth=0.8)
    if best_period:
        ax.axvline(x=best_period, color='red', linestyle='--', linewidth=1, label=f'Best P = {best_period:.4f} d')
        ax.legend()
    ax.set_xlabel("Period (days)")
    ax.set_ylabel("Power")
    ax.set_title(title)
    ax.minorticks_on()
    if save_path:
        fig.savefig(save_path)
        plt.close(fig)
    else:
        plt.show()
    return fig, ax
```

### 8. 天球投影（Mollweide）

```python
def plot_skymap(ra_deg, dec_deg, values=None, title="Sky Map", save_path=None):
    fig, ax = plt.subplots(figsize=(10, 5), subplot_kw={'projection': 'mollweide'})
    # 转换为弧度，RA 中心在 180°
    ra_rad = np.deg2rad(np.array(ra_deg))
    dec_rad = np.deg2rad(np.array(dec_deg))
    ra_rad = ra_rad - np.pi  # 中心化
    if values is not None:
        sc = ax.scatter(ra_rad, dec_rad, c=values, s=5, cmap='viridis', alpha=0.7, edgecolors='none')
        fig.colorbar(sc, ax=ax, pad=0.05, shrink=0.6)
    else:
        ax.scatter(ra_rad, dec_rad, c='black', s=5, alpha=0.5, edgecolors='none')
    ax.set_title(title)
    ax.grid(True, alpha=0.3)
    if save_path:
        fig.savefig(save_path, bbox_inches='tight')
        plt.close(fig)
    else:
        plt.show()
    return fig, ax
```

### 9. 发现图（Finding Chart）

```python
def plot_finding_chart(image_data, ra, dec, label=None, save_path=None):
    """image_data: 2D numpy array (from FITS or JPEG)"""
    fig, ax = plt.subplots()
    ax.imshow(image_data, cmap='gray_r', origin='lower')
    # 标注中心
    h, w = image_data.shape
    ax.plot(w/2, h/2, 'r+', markersize=15, markeredgewidth=2)
    if label:
        ax.annotate(label, (w/2, h/2), textcoords="offset points", xytext=(10, 10),
                    color='red', fontsize=10, fontweight='bold')
    ax.set_title(f"Finding Chart ({ra:.4f}, {dec:.4f})")
    ax.set_xlabel("X (pix)")
    ax.set_ylabel("Y (pix)")
    if save_path:
        fig.savefig(save_path)
        plt.close(fig)
    else:
        plt.show()
    return fig, ax
```

## 工作流程

1. **确定图表类型**：根据数据特征选择合适的图表模板
2. **准备数据**：从文献、数据库或分析结果中提取数据
3. **生成脚本**：使用上述模板，填入实际数据
4. **执行**：`python3 script.py`
5. **验证**：检查输出图片是否正确
6. **保存**：保存到 `figures/` 目录，命名格式：`{description}_{source}.{ext}`

## 保存路径约定

- `figures/` — 项目根目录下的图表输出目录
- 文件名格式：`{描述}_{来源}_{日期}.{格式}`
  - 示例：`spectrum_lamost_438809089_20260710.pdf`
  - 示例：`hr_diagram_m31_field_20260710.png`