<jupyter_start><jupyter_text>Emission line analysis In this tutorial we will analyze the Pa$\beta$ emission line in the $J$ band spectrum of GQ Lup B. This hydrogen emission line is related to the accretion processes in the circumsubstellar disk around this young, directly imaged brown dwarf.  Getting started We start by importing the required Python modules.<jupyter_code>import urllib.request
import species<jupyter_output><empty_output><jupyter_text>Next, we initiate the *species* workflow by running the [SpeciesInit](https://species.readthedocs.io/en/latest/species.core.html?highlight=speciesinitspecies.core.init.SpeciesInit) class.<jupyter_code>species.SpeciesInit()<jupyter_output>Initiating species v0.5.0... [DONE]
Creating species_config.ini... [DONE]
Database: /Users/tomasstolker/applications/species/docs/tutorials/species_database.hdf5
Data folder: /Users/tomasstolker/applications/species/docs/tutorials/data
Working folder: /Users/tomasstolker/applications/species/docs/tutorials
Creating species_database.hdf5... [DONE]
Creating data folder... [DONE]<jupyter_text>We will also download the VLT/SINFONI $J$ band spectrum of GQ Lup B that has been published by [Seifahrt et al. (2007)](https://ui.adsabs.harvard.edu/abs/2007A%26A...463..309S/abstract).<jupyter_code>urllib.request.urlretrieve('https://home.strw.leidenuniv.nl/~stolker/species/gqlupb_sinfoni_j.dat',
                           'gqlupb_sinfoni_j.dat')<jupyter_output><empty_output><jupyter_text>Adding data of GQ Lup B to the database To add data to the HDF5 database, we first create an instance of [Database](https://species.readthedocs.io/en/latest/species.data.html?highlight=databasespecies.data.database.Database). This object provides read and write access to the HDF5 file.<jupyter_code>database = species.Database()<jupyter_output><empty_output><jupyter_text>Next, we will add the available [magnitudes and distance](https://github.com/tomasstolker/species/blob/master/species/data/companions.py) of GQ Lup B to the database (only the distance is required later on for calculating the line luminosity) by using the [add_companion](https://species.readthedocs.io/en/latest/species.data.html?highlight=databasespecies.data.database.Database.add_companion) method of [Database](https://species.readthedocs.io/en/latest/species.data.html?highlight=databasespecies.data.database.Database). This method will add the magnitudes and filter profiles, but it also downloads a flux-calibrated spectrum of Vega and converts the magnitudes into fluxes.<jupyter_code>database.add_companion('GQ Lup B')<jupyter_output>Downloading Vega spectrum (270 kB)... [DONE]
Adding Vega spectrum... [DONE]
Adding filter: HST/WFPC2-PC.F606W... [DONE]
Adding filter: HST/WFPC2-PC.F814W... [DONE]
Adding filter: HST/NICMOS2.F171M... [DONE]
Adding filter: HST/NICMOS2.F190N... [DONE]
Adding filter: HST/NICMOS2.F215N... [DONE]
Adding filter: Magellan/VisAO.ip... [DONE]
Adding filter: Magellan/VisAO.zp... [DONE]
Adding filter: Magellan/VisAO.Ys... [DONE]
Adding filter: MKO/NSFCam.H... [DONE]
Adding filter: Paranal/NACO.NB405... [DONE]
Adding filter: Paranal/NACO.Mp... [DONE]
Adding filter: Paranal/NACO.Ks... [DONE]
Adding filter: Subaru/CIAO.CH4s... [DONE]
Adding filter: Subaru/CIAO.K... [DONE]
Adding filter: Subaru/CIAO.Lp... [DONE]
Adding object: GQ Lup B
   - Distance (pc) = 154.10 +/- 0.69
   - HST/WFPC2-PC.F606W:
      - Apparent magnitude = 19.19 +/- 0.07
      - Flux (W m-2 um-1) = 5.94e-16 +/- 3.83e-17
   - HST/WFPC2-PC.F814W:
      - Apparent magnitude = 17.67 +/- 0.05
      - Flux (W m-2 um-1) = 1.02e-15 +/- 4.6[...]<jupyter_text>Now that the photometric data of GQ Lup B are stored in the database, we will add the $J$ band spectrum by using the [add_object](https://species.readthedocs.io/en/latest/species.data.html?highlight=databasespecies.data.database.Database.add_object) method. Here, it is important that we use the same argument for `object_name` as was used as `name` with [add_companion](https://species.readthedocs.io/en/latest/species.data.html?highlight=databasespecies.data.database.Database.add_companion) (i.e. *GQ Lup B*). The argument of `spectrum` contains a dictionary with the spectra. Each value is a tuple with the spectrum, covariance matrix (set to `None` in this case), and the spectral resolution.<jupyter_code>database.add_object('GQ Lup B',
                    app_mag=None,
                    flux_density=None,
                    spectrum={'SINFONI': ('gqlupb_sinfoni_j.dat', None, 2500.)},
                    deredden=None)<jupyter_output>Adding object: GQ Lup B
   - Spectrum:
      - Database tag: SINFONI
      - Filename: gqlupb_sinfoni_j.dat
      - Data shape: (1014, 3)
      - Wavelength range (um): 1.20 - 1.35
      - Mean flux (W m-2 um-1): 3.27e-15
      - Mean error (W m-2 um-1): 2.79e-17
   - Spectral resolution:
      - SINFONI: 2500.0<jupyter_text>Initiating the line analysis We will now analyze the emission line by first creating an instance of [EmissionLine](https://species.readthedocs.io/en/latest/species.analysis.htmlspecies.analysis.emission_line.EmissionLine). The spectrum is read from the database by providing the same object name and spectrum name as were used with [add_object](https://species.readthedocs.io/en/latest/species.data.html?highlight=databasespecies.data.database.Database.add_object). We will select a limited wavelength range around the Pa$\beta$ line by setting the argument of `wavel_range`.<jupyter_code>line = species.EmissionLine(object_name='GQ Lup B',
                            spec_name='SINFONI',
                            lambda_rest=1.282,
                            wavel_range=(1.26, 1.32))<jupyter_output><empty_output><jupyter_text>Let's have a look at the data of the cropped $J$ spectrum that is extracted.<jupyter_code>line.spectrum<jupyter_output><empty_output><jupyter_text>Subtracting the continuum flux In the $J$ band, we need to subtract the pseudo-continuum of the atmospheric emission. To do so, we make use of the [subtract_continuum](https://species.readthedocs.io/en/latest/species.analysis.htmlspecies.analysis.emission_line.EmissionLine.subtract_continuum) method of [EmissionLine](https://species.readthedocs.io/en/latest/species.analysis.htmlspecies.analysis.emission_line.EmissionLine). The continuum is modeled with a 3rd order polynomial (after smoothing the spectrum) and we use a linear least square fitting algorithm to optimize the parameters. A plot with the best-fit polynomial and the continuum-subtracted spectrum is created, which can be inspected before continuing with the line analysis.<jupyter_code>line.subtract_continuum(poly_degree=3,
                        plot_filename=None)<jupyter_output>Fitting continuum... [DONE]
Plotting continuum fit... [DONE]<jupyter_text>Integrating the line flux Now that the continuum is subtracted, we can determine the line characteristics in two different ways. Here, we will measure the line properties directly from the spectrum. Later on, we will fit a Gaussian profile and adopt the sampled parameters for inferring the line properties.We use the [integrate_flux](https://species.readthedocs.io/en/latest/species.analysis.htmlspecies.analysis.emission_line.EmissionLine.integrate_flux) method to integrate the line flux across the specified wavelength range of ``wavel_int`` after interpolating the spectrum to $R = 100000$. Several other characteristics are also calculated from the spectral line. The uncertainties are estimated by sampling 1000 random spectra from the errors on the fluxes.<jupyter_code>line.integrate_flux(wavel_int=(1.280, 1.285),
                    interp_kind='linear',
                    plot_filename=None)<jupyter_output>Plotting integrated line... [DONE]<jupyter_text>Modeling the Pa$\beta$ emission line We will now use the [fit_gaussian](https://species.readthedocs.io/en/latest/species.analysis.htmlspecies.analysis.emission_line.EmissionLine.fit_gaussian) method of [EmissionLine](https://species.readthedocs.io/en/latest/species.analysis.htmlspecies.analysis.emission_line.EmissionLine) to fit the spectrum with a Gaussian function while mapping the posterior distributions of the 3 free parameters (amplitude, mean, and standard deviation) with the nested sampling algorithm of [UltraNest](https://johannesbuchner.github.io/UltraNest/index.html). The argument of `bounds` is a dictionary with the prior boundaries of the three parameters (`gauss_amplitude`, `gauss_mean`, and `gauss_sigma`). If any parameter is missing or `bounds=None` then conservative boundaries will be estimated from the spectrum. The posterior samples will be stored in the database with the `tag` name and a plot is created with a comparison of the data and best-fit model (i.e. the median parameter values).<jupyter_code>line.fit_gaussian(tag='p-beta',
                  min_num_live_points=400,
                  bounds={'gauss_amplitude': (0., 2e-15)},
                  output='ultranest',
                  plot_filename=None,
                  show_status=False)<jupyter_output>Creating directory for new run ultranest/run1
[ultranest] Sampling 400 live points from prior ...<jupyter_text>Plotting the posterior distributions Now that the posterior samples are available in the database, we can use the [plot_posterior](https://species.readthedocs.io/en/latest/species.plot.html?highlight=plot_posteriorspecies.plot.plot_mcmc.plot_posterior) function for visualizing the 1D and 2D projected distributions with [corner.py](https://corner.readthedocs.io). The argument of `tag` is set to the same database tag that was used with [fit_gaussian](https://species.readthedocs.io/en/latest/species.analysis.htmlspecies.analysis.emission_line.EmissionLine.fit_gaussian).<jupyter_code>species.plot_posterior(tag='p-beta',
                       burnin=None,
                       title=None,
                       title_fmt=['.2e', '.3f', '.3f', '.1f', '.2e', '.2e', '.2f', '.0f'],
                       offset=(-0.4, -0.35),
                       output=None)<jupyter_output>Median sample:
   - gauss_amplitude = 1.23e-15
   - gauss_mean = 1.28e+00
   - gauss_sigma = 4.29e-04
   - gauss_fwhm = 2.36e+02
   - line_flux = 1.32e-18
   - line_luminosity = 9.81e-07
   - line_eq_width = -3.72e+00
   - line_vrad = 1.50e+01
   - distance = 1.54e+02
Plotting the posterior... [DONE]
StarcoderdataJupyterStructuredCleanDedup