The incremental functions split up the spectrum into small sections an calculate the precision on each. This enables the change in each precision/quality accross the whole wavelength range to be observed. These were inspired by the plots in Bouchy 2001 and Artigau 2018.
Here we show how to use easily use the incremental_quality and incremental_rv functions.
The split the specturm in sections with a width of X% the current wavelength. This is logarithmic and the section size increases with wavelength.
In [1]:
import matplotlib.pyplot as plt
from eniric.precision import incremental_rv, incremental_quality
from scripts.phoenix_precision import convolve_and_resample
from eniric.snr_normalization import snr_constant_band
from eniric.utilities import load_aces_spectrum
In [2]:
# Settings
band = "J" # analyse this band only (for speed)
snr = 100 # SNR per pixel normalization at center of band.
R = 100000 # Resolution
vsini = 1 # Rotation speed per pixel
smpl=3 # Sampling
In [3]:
wav1, flux1 = load_aces_spectrum([3900, 4.5, 0.0, 0])
wav2, flux2 = load_aces_spectrum([3500, 4.5, 0.0, 0])
wav3, flux3 = load_aces_spectrum([2800, 4.5, 0.0, 0])
wav4, flux4 = load_aces_spectrum([2600, 4.5, 0.0, 0])
In [4]:
# Convolution and resampling
wav1, flux1 = convolve_and_resample(
wav1, flux1, vsini=vsini, R=R, band=band, sampling=smpl
)
wav2, flux2 = convolve_and_resample(
wav2, flux2, vsini=vsini, R=R, band=band, sampling=smpl
)
wav3, flux3 = convolve_and_resample(
wav3, flux3, vsini=vsini, R=R, band=band, sampling=smpl
)
wav4, flux4 = convolve_and_resample(
wav4, flux4, vsini=vsini, R=R, band=band, sampling=smpl
)
In [5]:
# The 4 spectra examples
plt.plot(wav1, flux1, label="M0")
plt.plot(wav2, flux2, label="M3")
plt.plot(wav3, flux3, label="M6")
plt.plot(wav4, flux4, label="M9")
plt.legend()
plt.ylabel("Flux")
plt.xlabel("Wavelength ($\mu$m)")
plt.show()
In [6]:
# Only dealing with one band so band=band
norm1 = snr_constant_band(wav1, flux1, band=band, snr=snr)
norm2 = snr_constant_band(wav2, flux2, band=band, snr=snr)
norm3 = snr_constant_band(wav3, flux3, band=band, snr=snr)
norm4 = snr_constant_band(wav4, flux4, band=band, snr=snr)
flux1 /= norm1
flux2 /= norm2
flux3 /= norm3
flux4 /= norm4
In [7]:
print("This should now be 100! if band == J")
norm4 = snr_constant_band(wav4, flux4, band=band, snr=snr)
In [8]:
x1, q1 = incremental_quality(wav1, flux1, percent=1)
x2, q2 = incremental_quality(wav2, flux2, percent=1)
x3, q3 = incremental_quality(wav3, flux3, percent=1)
x4, q4 = incremental_quality(wav4, flux4, percent=1)
In [9]:
plt.plot(x1, q1, label="M0", drawstyle="steps-mid")
plt.plot(x2, q2, label="M3", drawstyle="steps-mid")
plt.plot(x3, q3, label="M6", drawstyle="steps-mid")
plt.plot(x4, q4, label="M9", drawstyle="steps-mid")
plt.legend()
plt.ylabel("Spectral Quality")
plt.xlabel("Wavelength ($\mu$m)")
plt.show()
In [10]:
xrv1, rv1 = incremental_rv(wav1, flux1, mask=None, percent=1)
xrv2, rv2 = incremental_rv(wav2, flux2, mask=None, percent=1)
xrv3, rv3 = incremental_rv(wav3, flux3, mask=None, percent=1)
xrv4, rv4 = incremental_rv(wav4, flux4, mask=None, percent=1)
In [11]:
plt.plot(xrv1, rv1, label="M0", drawstyle="steps-mid")
plt.plot(xrv2, rv2, label="M3", drawstyle="steps-mid")
plt.plot(xrv3, rv3, label="M6", drawstyle="steps-mid")
plt.plot(xrv4, rv4, label="M9", drawstyle="steps-mid")
plt.legend()
plt.ylabel("RV precision (m/s)")
plt.xlabel("Wavelength ($\mu$m)")
plt.show()
In [ ]: