compare_model_spectra_ASPCAP.ipynb Enhancement

https://github.com/sdss/dr20_tutorials/blob/main/compare_model_spectra_ASPCAP.ipynb

Add enhancement code:

------------------------------------------------------------

NEW CELL: Auto-detect each element's diagnostic window + zoomed residual plot

------------------------------------------------------------

'''
Addresses the notebook's own stated limitation (Cell 16): the full-width
plots hide line-level detail. Rather than guessing exact line-center
wavelengths from memory, this finds where each element's model_flux
differs most from the full best-fit model -- those pixels are, by
construction, the ASPCAP fitting window for that element -- and zooms
into a region around them, with a residual panel underneath.
'''
full_model = aspcapStar["model_flux"][0]
window_half_width = 15 # Angstroms of context around the detected window

for elem in elems:
elem_model = aspcapStar[f"model_flux_{elem}_h"][0]
diff = np.abs(elem_model - full_model)

if not np.any(np.isfinite(diff)) or np.nanmax(diff) == 0:
    print(f"{elem}: no detectable difference from full model, skipping")
    continue

# locate the strongest single diagnostic feature for this element
peak_idx = np.nanargmax(diff)
center_wl = wavelength[peak_idx]
window = (wavelength > center_wl - window_half_width) & (wavelength < center_wl + window_half_width)
obs_window = (obs_wl > center_wl - window_half_width) & (obs_wl < center_wl + window_half_width) & mask

fig, (ax1, ax2) = plt.subplots(nrows=2, figsize=(9, 6), sharex=True,
                                gridspec_kw={'height_ratios': [3, 1]})
ax1.plot(obs_wl[obs_window], obs_flux[obs_window], color='k', label='observed', lw=1.2)
ax1.plot(wavelength[window], (full_model * continuum)[window], '--', label='full best-fit model')
ax1.plot(wavelength[window], (elem_model * continuum)[window], '--',
          label=f"{elem.title()} model")
ax1.axvline(center_wl, color='gray', lw=0.5, ls=':')
ax1.set_ylabel("flux")
ax1.set_title(fr"{elem.title()}: [{elem.title()}/H]={aspcapStar[f'{elem}_h'][0]:.3f} "
              fr"$\pm$ {aspcapStar[f'e_{elem}_h'][0]:.3f}, zoomed on strongest diagnostic pixel")
ax1.legend(fontsize=8)

# residual: element model vs full model, in normalized (continuum) units
ax2.plot(wavelength[window], (elem_model - full_model)[window], color='firebrick')
ax2.axhline(0, color='k', lw=0.7, ls='--')
ax2.set_ylabel("elem - full\nmodel")
ax2.set_xlabel("λ [Å]")
plt.tight_layout()
plt.show()