{ "cells": [ { "cell_type": "markdown", "id": "8c8138ef-e687-444e-95b1-9797654f3f62", "metadata": {}, "source": [ "# Spectral Axes: Velocity Frames and Doppler Conventions\n", "\n", "In this notebook we show how to change the spectral axis to represent different rest frames and Doppler conventions. In addition, a spectrum is normally\n", "associated with a redshift/radial velocity, which can also be modified. We also show how to save a spectrum in\n", "frequency and velocity space.\n", "\n", "Supported **doppler conventions**: (see also https://www.gb.nrao.edu/~fghigo/gbtdoc/doppler.html)\n", "1. radio\n", "2. relativistic\n", "3. optical\n", "\n", "Some of the supported **rest frames**:\n", "1. itrs - topocentric\n", "2. icrs - barycentric\n", "3. gcrs - geocentric\n", "4. hcrs - heliocentric\n", "5. lsrk, lsrd - Local Standard of Rest (Kinematic or Dynamic)\n", "\n", "Although the `Spectrum.plot()` can plot a spectrum in different frames with different conventions, it does not modify the underlying data and meta-data in the spectrum. To make these persistent (e.g. necessary when writing a spectrum) there is both an in-place and copy operation to modify frame and convention. Here's a dysh command summary that we will cover in this notebook, leaving out the values and arguments that are not relevant:\n", "\n", "```\n", " ta.plot(vel_frame=, doppler_convention=)\n", "\n", " ta.velocity_axis_to(toframe=, doppler_convention=)\n", "\n", " ta.set_frame()\n", " ta.set_convention()\n", "\n", " ta1 = ta.with_frame()\n", " ta2 = ta1.with_velocity_convention()\n", "\n", " ta3 = ta.with_spectral_axis_unit(\"km/s\")\n", "\n", " ta.set_redshift.to()\n", " ta.set_radial_velocity_to()\n", " ta.shift_spectrum_to(redshift=)\n", " ta.shift_spectrum_to(radial_velocity=)\n", "\n", " ta.rest_value =\n", "\n", " # no setters for these\n", " ta.redshift\n", " ta.radial_velocity\n", "```\n", "\n", "\n", "## Loading Modules\n", "We start by loading the modules we will use for the data reduction. \n" ] }, { "cell_type": "code", "execution_count": 1, "id": "b4967550-2ca1-4931-b53b-6f9868718490", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "import" ] }, "outputs": [], "source": [ "# These modules are required for working with the data.\n", "from dysh.fits.gbtfitsload import GBTFITSLoad\n", "from dysh.log import init_logging\n", "from astropy import units as u\n", "from dysh.spectra.spectrum import Spectrum\n", "\n", "# These modules are used for file I/O\n", "from dysh.util.files import dysh_data\n", "from pathlib import Path" ] }, { "cell_type": "markdown", "id": "c91993b7-fc23-4164-82a6-d8446cf7f5af", "metadata": {}, "source": [ "## Setup\n", "We start the dysh logging, so we get more information about what is happening.\n", "This is only needed if working on a notebook.\n", "If using the CLI through the dysh command, then logging is setup for you." ] }, { "cell_type": "code", "execution_count": 2, "id": "4a3ea33e-0d69-49ba-9aaf-74ac4815f1c7", "metadata": {}, "outputs": [], "source": [ "init_logging(2)\n", "\n", "# also create a local \"output\" directory where temporary notebook files can be stored.\n", "output_dir = Path.cwd() / \"output\"\n", "output_dir.mkdir(exist_ok=True)" ] }, { "cell_type": "markdown", "id": "87669763-8d96-4521-9e81-f69d7213e133", "metadata": {}, "source": [ "## Data Retrieval\n", "\n", "Download the example SDFITS data, if necessary." ] }, { "cell_type": "code", "execution_count": 3, "id": "6bc88bc5-986d-4eae-b1c7-6398cc9ddd5a", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "wget" ] }, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "20:22:11.259 I Resolving test=getps -> AGBT05B_047_01/AGBT05B_047_01.raw.acs/\n" ] } ], "source": [ "filename = dysh_data(test=\"getps\")" ] }, { "cell_type": "markdown", "id": "e05facfd-26bc-484c-8d35-a2f3fa73086d", "metadata": {}, "source": [ "## Data Loading" ] }, { "cell_type": "code", "execution_count": 4, "id": "93a62e3a-c95d-475b-8602-b5b8b7934733", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "load" ] }, "outputs": [ { "data": { "text/html": [ "\n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", " \n", "
SCANOBJECTVELOCITYPROCPROCSEQNRESTFREQDOPFREQ# IF# POL# INT# FEEDAZIMUTHELEVATION
51NGC52914386.0OnOff11.4204051.42040512111198.343118.6427
52NGC52914386.0OnOff21.4204051.42040512111198.930618.7872
53NGC52914386.0OnOff11.4204051.42040512111199.330518.3561
54NGC52914386.0OnOff21.4204051.42040512111199.915718.4927
55NGC52914386.0OnOff11.4204051.42040512111200.304218.0575
56NGC52914386.0OnOff21.4204051.42040512111200.890618.1860
57NGC52914386.0OnOff11.4204051.42040512111202.327517.3853
58NGC52914386.0OnOff21.4204051.42040512111202.919217.4949
" ], "text/plain": [ "" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "sdfits = GBTFITSLoad(filename)\n", "sdfits.summary()" ] }, { "cell_type": "markdown", "id": "d3ef13a3-0037-4a75-8011-bb01cb9cb524", "metadata": {}, "source": [ "More background on the GBT SDFITS files can be found on\n", "https://dysh.readthedocs.io/en/latest/reference/sdfits_files/gbt_sdfits.html" ] }, { "cell_type": "markdown", "id": "e84c5009-96f4-489d-bc9d-89de390d53c9", "metadata": {}, "source": [ "## Data Reduction\n", "\n", "Next we fetch and calibrate the position switched data. We will use this data to show how to change rest frames and Doppler conventions. All of these occur at the spectrum level, so we need to get a spectrum first. \n", "More details can be found in the \n", "[Position Switching](https://dysh.readthedocs.io/en/latest/users_guide/positionswitch.html)\n", "notebook.\n", "\n", "We use the time-averaged spectrum from scans 51/52:" ] }, { "cell_type": "code", "execution_count": 5, "id": "e0a33c19-576b-46e8-9a57-f34b8f85f6aa", "metadata": { "editable": true, "slideshow": { "slide_type": "" }, "tags": [ "timeaverage", "test" ] }, "outputs": [], "source": [ "ta = sdfits.getps(scan=51, ifnum=0, plnum=0, fdnum=0).timeaverage()" ] }, { "cell_type": "markdown", "id": "d8959e6e-90f4-4af5-a936-1c36af7c0186", "metadata": {}, "source": [ "## Changing the x-axis of a `Spectrum` Plot\n", "Note this changes the axis of the plot but does not affect the underlying Spectrum object.\n" ] }, { "cell_type": "markdown", "id": "94f52f21-0af6-4e6f-9be5-50d37e1f5563", "metadata": {}, "source": [ "### Default Rest Frame\n", "\n", "The default plot uses the frequency frame and Doppler convention found in the SDFITS file. In this case, that is topocentric frame (ITRS) and the optical convention. Notice that the GBT SDFITS file stores frequencies with CTYPE1=\"FREQ-OBS\". Check by looking at sdfits[\"CTYPE1\"][0]." ] }, { "cell_type": "code", "execution_count": 6, "id": "dd91a385-54f3-4cde-8e90-66f06dc93ed4", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "'FREQ-OBS'" ] }, "execution_count": 6, "metadata": {}, "output_type": "execute_result" } ], "source": [ "sdfits[\"CTYPE1\"][0]" ] }, { "cell_type": "code", "execution_count": 7, "id": "a7b7e8ec-0ef0-45e1-817a-e4aebec4330a", "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "f854b3e35cb141efbb02d93a0902ab71", "version_major": 2, "version_minor": 0 }, "text/plain": [ "VBox(children=(HBox(children=(Button(description='Clear All Regions', style=ButtonStyle(), tooltip='Clear all …" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "ta.plot();" ] }, { "cell_type": "markdown", "id": "adbb99c0-f0d5-430b-a4a2-11969102ccb9", "metadata": {}, "source": [ "Although the spectrum `ta` can be printed using `print`, just printing it as an object gives some more information" ] }, { "cell_type": "code", "execution_count": 8, "id": "c5dc53c4-c1b3-4a36-b66e-aeffae6b872c", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Spectrum (length=32768)\n", "Flux=[0.2448421 0.31819268 0.19866335 ... 0.57650381 0.2179878\n", " 1.20767879] K, mean=0.28071 K\n", "Spectral Axis=[1.42481684e+09 1.42481531e+09 1.42481379e+09 ...\n", " 1.37482142e+09 1.37481989e+09 1.37481836e+09] Hz, mean=1399817601.06055 Hz\n", "\n", "\n" ] }, { "data": { "text/plain": [ "\n", " target: \n", " observer to target (computed from above):\n", " radial_velocity=4410.070039086893 km / s\n", " redshift=0.014820217767934851\n", " doppler_rest=1420405000.0 Hz\n", " doppler_convention=optical)\n", " [1.42481684e+09 1.42481531e+09 1.42481379e+09 ... 1.37482142e+09\n", " 1.37481989e+09 1.37481836e+09] Hz> (length=32768))>" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "print(ta)\n", "print(\"\\n\")\n", "ta" ] }, { "cell_type": "markdown", "id": "d9834355-10fc-40f8-9099-5d51d267a07c", "metadata": {}, "source": [ "A large amount of information is also stores in the `meta` data, which is a python dictionary associated with the Spectrum.\n", "Here is how to access the meta data of a spectrum, where we use a python trick to make the dictionary shown in alphabetical order:" ] }, { "cell_type": "code", "execution_count": 9, "id": "6353c82c-9711-4580-85d9-9faf1191439d", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "{'AP_EFF': np.float64(0.7047451620556072),\n", " 'AZIMUTH': 198.1588901992111,\n", " 'BACKEND': 'Spectrometer',\n", " 'BANDWID': 50000000.0,\n", " 'BINTABLE': 0,\n", " 'BUNIT': 'K',\n", " 'CAL': 'F',\n", " 'CALPOSITION': 'Unknown',\n", " 'CALTYPE': 'LOW',\n", " 'CDELT1': -1525.87890625,\n", " 'CRPIX1': 16385.0,\n", " 'CRVAL1': 1399816838.1210938,\n", " 'CRVAL2': np.float64(206.85210757719534),\n", " 'CRVAL3': np.float64(-30.407015310885345),\n", " 'CRVAL4': -6,\n", " 'CTYPE1': 'FREQ-OBS',\n", " 'CTYPE2': 'RA',\n", " 'CTYPE3': 'DEC',\n", " 'CTYPE4': 'STOKES',\n", " 'CUNIT1': 'Hz',\n", " 'CUNIT2': 'deg',\n", " 'CUNIT3': 'deg',\n", " 'DATE': '2024-03-15T14:05:45',\n", " 'DATE-OBS': '2005-06-27T02:05:58.00',\n", " 'DOPFREQ': 1420405000.0,\n", " 'DURATION': np.float64(55.5225),\n", " 'ELEVATIO': 18.69449298449073,\n", " 'EQUINOX': 2000.0,\n", " 'EXPOSURE': np.float64(53.7157844463702),\n", " 'EXTEND': True,\n", " 'EXTNAME': 'SINGLE DISH',\n", " 'FDNUM': 0,\n", " 'FEED': 1,\n", " 'FEEDEOFF': 0.0,\n", " 'FEEDXOFF': 0.0,\n", " 'FITSINDEX': 0,\n", " 'FITSVER': '1.9',\n", " 'FREQRES': 1846.3134765625,\n", " 'FRONTEND': 'Rcvr1_2',\n", " 'HDU': 1,\n", " 'HUMIDITY': 0.754,\n", " 'IFNUM': 0,\n", " 'INSTRUME': 'Spectrometer',\n", " 'LASTOFF': 0,\n", " 'LASTON': 51,\n", " 'LST': 54470.49953792763,\n", " 'MEANTSYS': np.float64(19.356765383560564),\n", " 'NAXIS1': 32768,\n", " 'OBJECT': 'NGC5291',\n", " 'OBSERVER': 'Jeff Mangum',\n", " 'OBSFREQ': 1399816838.1210938,\n", " 'OBSID': 'unknown',\n", " 'OBSMODE': 'OnOff:PSWITCHON:TPWCAL',\n", " 'OBSTYPE': 'PSWITCHON',\n", " 'ORIGIN': 'NRAO Green Bank',\n", " 'PLNUM': 0,\n", " 'PRESSURE': 697.7295162882529,\n", " 'PROC': 'OnOff',\n", " 'PROCSCAN': 'Unknown',\n", " 'PROCSEQN': 1,\n", " 'PROCSIZE': 2,\n", " 'PROCTYPE': 'SIMPLE',\n", " 'PROJID': 'AGBT05B_047_01',\n", " 'QD_BAD': -1,\n", " 'QD_EL': nan,\n", " 'QD_METHOD': '',\n", " 'QD_XEL': nan,\n", " 'RADESYS': 'FK5',\n", " 'RESTFREQ': 1420405000.0,\n", " 'RESTFRQ': np.float64(1420405000.0),\n", " 'ROW': 3,\n", " 'RVSYS': 4376523.966595529,\n", " 'SAMPLER': 'A13',\n", " 'SCAN': 51,\n", " 'SDFITVER': 'sdfits ver1.22',\n", " 'SE_UNIT': 'micron',\n", " 'SIDEBAND': 'U',\n", " 'SIG': 'T',\n", " 'SIMPLE': True,\n", " 'SITEELEV': 824.595,\n", " 'SITELAT': 38.43312,\n", " 'SITELONG': -79.83983,\n", " 'SRFEED': 0,\n", " 'SUBOBSMODE': 'TPWCAL',\n", " 'SUBREF_STATE': 1,\n", " 'SURF_ERR': np.float64(229.99999999999997),\n", " 'TAMBIENT': 292.54999999999995,\n", " 'TCAL': np.float64(1.4526499509811401),\n", " 'TCOLD': nan,\n", " 'TDIM7': '(32768,1,1,1)',\n", " 'TELESCOP': 'NRAO_GBT',\n", " 'TIMESTAMP': '2005_06_27_02:05:58',\n", " 'TRGTLAT': -30.407,\n", " 'TRGTLONG': 206.85199999999995,\n", " 'TSCALE': 'Ta',\n", " 'TSCALFAC': np.float64(1.0),\n", " 'TSYS': np.float64(19.353858908450768),\n", " 'TUNIT7': 'K',\n", " 'TWARM': nan,\n", " 'VELDEF': 'OPTI-LSR',\n", " 'VELOCITY': 4386000.0,\n", " 'VFRAME': 22609.232181632025,\n", " 'WTTSYS': np.float64(19.353858908450764),\n", " 'ZEROCHAN': 282.5385437011719}" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "dict(sorted(ta.meta.items()))" ] }, { "cell_type": "markdown", "id": "4d97f4ae-3fbd-4561-94eb-df3f07d27ff5", "metadata": {}, "source": [ "For the remainder of the notebook, we focus on the small section of the spectrum near 1.39 GHz where a random spike can be seen in the previous plot. After some trial and error this is channel 21921 where the peak value is 0.59789101 K.\n", "\n", "Let's plot this zoomed spectrum in km/s so we can better see the effect of changing frames and conventions." ] }, { "cell_type": "code", "execution_count": 10, "id": "e894fb31-87a2-4204-8268-83bf12520e04", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Taking a zoomed spectrum\n", "Spectrum (length=3)\n", "Flux=[0.33822484 0.59789101 0.26609857] K, mean=0.40074 K\n", "Spectral Axis=[1.39136957e+09 1.39136805e+09 1.39136652e+09] Hz, mean=1391368046.61719 Hz\n" ] }, { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "9b8e4370269b45e99fe95f553c642f9b", "version_major": 2, "version_minor": 0 }, "text/plain": [ "VBox(children=(HBox(children=(Button(description='Clear All Regions', style=ButtonStyle(), tooltip='Clear all …" ] }, "metadata": {}, "output_type": "display_data" }, { "name": "stdout", "output_type": "stream", "text": [ "Peak at 6256.475163870069 km / s in the itrs reference frame\n" ] } ], "source": [ "# pick a zoomed spectrum (this cell is safe to re-execute)\n", "if ta.nchan > 1000:\n", " print(f\"Taking a zoomed spectrum\")\n", " print(ta[21920:21923]) # print values around the spike, peak is at channel 21921\n", " ta = ta[21850:22000] # these are the 150 channels we will zoom into\n", "else:\n", " print(f\"Already had the zoomed spectrum of {ta.nchan} channels\")\n", " print(ta[70:73]) # peak is at channel 71\n", "\n", "ta.plot(xaxis_unit='km/s')\n", "print(f\"Peak at {ta.velocity_axis_to('km/s')[71]} in the {ta.velocity_frame} reference frame\")" ] }, { "cell_type": "markdown", "id": "33276c88-f6a1-48a2-bb0c-6a14803b6daf", "metadata": {}, "source": [ "## Writing a spectrum \n", "\n", "The default \n", "[Spectrum.write](https://specutils.readthedocs.io/en/stable/api/specutils.Spectrum.html#specutils.Spectrum.write)\n", "uses its native units (Hz) to write a spectrum. Arguably more useful is the spectral axis in km/s. One way is to make a copy of the spectrum while changing the units, using \n", "[Spectrum.with_spectral_axis_unit](https://specutils.readthedocs.io/en/stable/api/specutils.Spectrum.html#specutils.Spectrum.with_spectral_axis_unit)\n", "\n", "But first a default write operation, with spectral axis in Hz:" ] }, { "cell_type": "code", "execution_count": 11, "id": "f2eb5da4-4219-44db-923d-f7d39b270a6a", "metadata": {}, "outputs": [], "source": [ "# helper function to show the first few lines of a file\n", "def head(filename, nlines=3):\n", " \"\"\" emulate the unix `head` program\n", " \"\"\"\n", " with open(filename, 'r') as file:\n", " for i, line in enumerate(file, 1):\n", " print(line.strip())\n", " if i == nlines:\n", " break\n", "\n", "spike_table = \"output/ngc5291_spike.tab\"" ] }, { "cell_type": "code", "execution_count": 12, "id": "8d6916ce-9e4c-433f-bb97-90c06091e4fd", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "# spectral_axis flux uncertainty weight mask baseline\n", "1391476384.0195374 0.34372934266958427 0.0 218.81991773344748 0 0.0\n", "1391474858.1406312 0.27963123045821414 0.0 218.81991773344748 0 0.0\n" ] } ], "source": [ "ta.write(spike_table,format=\"ascii.commented_header\", overwrite=True) \n", "head(spike_table)" ] }, { "cell_type": "code", "execution_count": 13, "id": "c7c383d5-1c2f-4838-b7db-d2531ab669a5", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "# spectral_axis flux uncertainty weight mask baseline\n", "6232.6468425347 0.34372934266958427 0.0 218.75420976172921 0 0.0\n", "6232.982426932423 0.27963123045821414 0.0 218.75420976172921 0 0.0\n" ] } ], "source": [ "ta1 = ta.with_spectral_axis_unit(\"km/s\")\n", "ta1.write(spike_table,format=\"ascii.commented_header\", overwrite=True) \n", "head(spike_table)" ] }, { "cell_type": "markdown", "id": "d83e0206-a79f-44c7-bd2b-a45516720e75", "metadata": {}, "source": [ "### Change Rest Frame\n", "You can change the velocity frame by supplying one of the [built-in astropy coordinate frames](https://docs.astropy.org/en/stable/coordinates/index.html#built-in-frame-classes). These are specified by a string name.\n", "\n", "For example, to plot in the barycentric frame use `\"icrs\"`. The change from topocentric to barycentric is small, about 25 km/s in this case." ] }, { "cell_type": "code", "execution_count": 14, "id": "11ccdc23-44c8-4bc7-b70b-9432166e2a66", "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "01517a3bd0a54f91a98f80ae81fd064c", "version_major": 2, "version_minor": 0 }, "text/plain": [ "VBox(children=(HBox(children=(Button(description='Clear All Regions', style=ButtonStyle(), tooltip='Clear all …" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "ta.plot(vel_frame='icrs', xaxis_unit='km/s', grid=True);" ] }, { "cell_type": "markdown", "id": "ae93d848-3b37-4853-a0eb-86733362d2d3", "metadata": {}, "source": [ "Recall that the peak was at 6256.5 km/s in the topocentric frame, whereas in the barycentric frame it is more like 6231.9 km/s." ] }, { "cell_type": "markdown", "id": "011f3b4f-2749-4baf-b25e-0d55c4943c29", "metadata": {}, "source": [ "In addition to the `astropy` frame names, we also allow `'topo'` and `'topocentric'`." ] }, { "cell_type": "code", "execution_count": 15, "id": "a63ccae9-0c5e-4d6a-bf6b-707e37663fef", "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "3b2ac03f54ae408cb04ebd071617e5db", "version_major": 2, "version_minor": 0 }, "text/plain": [ "VBox(children=(HBox(children=(Button(description='Clear All Regions', style=ButtonStyle(), tooltip='Clear all …" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "ta.plot(vel_frame='topo', xaxis_unit='km/s');" ] }, { "cell_type": "markdown", "id": "bda0af0e-35b1-4061-a2e1-b12827e3634c", "metadata": {}, "source": [ "### Doppler Convention \n", "One can also change the Doppler convention between `radio`, `optical`, and `relativistic`.\n", "Here we also change the x-axis to velocity units and use LSRK frame." ] }, { "cell_type": "code", "execution_count": 16, "id": "e66e698f-42dd-490a-9727-f6493be31670", "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "3acc426fe20247ba91f21446a9ce9b75", "version_major": 2, "version_minor": 0 }, "text/plain": [ "VBox(children=(HBox(children=(Button(description='Clear All Regions', style=ButtonStyle(), tooltip='Clear all …" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "ta.plot(vel_frame='lsrk', doppler_convention='radio', xaxis_unit='km/s');" ] }, { "cell_type": "markdown", "id": "73703ac1-dc37-4da6-b4bb-cacfe8ab2578", "metadata": {}, "source": [ "## Changing the Spectral Axis of the `Spectrum`.\n", "There are two ways to accomplish this. One returns a copy of the original Spectrum with the new spectral axis; \n", "the other changes the spectral axis in place.\n", "\n", "The `with_frame()` and `with_velocity_convention()` are used when returning a copy. The two can be chained.\n", "\n", "The `set_frame()` and `set_convention()` are used in place." ] }, { "cell_type": "markdown", "id": "4ad12c65-a95e-498f-98f5-6646356a0285", "metadata": {}, "source": [ "### A. Return a copy of the spectrum using [Spectrum.with_frame](https://dysh.readthedocs.io/en/latest/reference/modules/dysh.spectra.html#dysh.spectra.spectrum.Spectrum.with_frame)\n" ] }, { "cell_type": "code", "execution_count": 17, "id": "df14fc13-c5f0-46c5-9d4f-2c839c79a197", "metadata": { "scrolled": true }, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "f68d431e46a743a9b6d338d7ced8088e", "version_major": 2, "version_minor": 0 }, "text/plain": [ "VBox(children=(HBox(children=(Button(description='Clear All Regions', style=ButtonStyle(), tooltip='Clear all …" ] }, "metadata": {}, "output_type": "display_data" }, { "data": { "text/plain": [ "" ] }, "execution_count": 17, "metadata": {}, "output_type": "execute_result" } ], "source": [ "newspec = ta.with_frame('galactocentric')\n", "newspec.plot(xaxis_unit=\"km/s\",doppler_convention='radio')" ] }, { "cell_type": "code", "execution_count": 18, "id": "56c65a59-f7ed-4b1c-8086-097cf42878ea", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "The new spectral axis frame is galactocentric\n", "radio: 5976.071520424437 km / s\n", "default: 6097.621687333563 km / s optical\n" ] } ], "source": [ "print(f\"The new spectral axis frame is {newspec.velocity_frame}\")\n", "print('radio: ',newspec.velocity_axis_to(\"km/s\", doppler_convention=\"radio\")[71]) # 5976.1\n", "print('default:',newspec.velocity_axis_to(\"km/s\")[71], newspec.velocity_convention) # 6097.6" ] }, { "cell_type": "markdown", "id": "f9cd5c17-5e5f-412c-9f5e-d9a63f87c79d", "metadata": {}, "source": [ "One can see that the spectral axis of the new spectrum is different. About 722 kHz." ] }, { "cell_type": "code", "execution_count": 19, "id": "765c5f26-451f-40cf-8b03-d105a44d884a", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$[-722615.33,~-722614.54,~-722613.75,~\\dots,~-722498.85,~-722498.05,~-722497.26] \\; \\mathrm{Hz}$" ], "text/plain": [ "" ] }, "execution_count": 19, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ta.spectral_axis - newspec.spectral_axis\n" ] }, { "cell_type": "markdown", "id": "2f84bba9-13e2-434a-94ce-7bc1c66e0cd5", "metadata": {}, "source": [ "### B. Change the spectral axis in place using [Spectrum.set_frame](https://dysh.readthedocs.io/en/latest/reference/modules/dysh.spectra.html#dysh.spectra.spectrum.Spectrum.set_frame)\n" ] }, { "cell_type": "code", "execution_count": 20, "id": "11b429f1-89c7-47e2-91e3-9057e22d5e7b", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "6256.475163870069 km / s\n", "Changed spectral axis frame to gcrs\n", "6256.365280719197 km / s\n" ] } ], "source": [ "print(ta.velocity_axis_to(\"km/s\")[71]) # 6256.475\n", "sa = ta.spectral_axis\n", "ta.set_frame('gcrs')\n", "print(f\"Changed spectral axis frame to {ta.velocity_frame}\")\n", "print(ta.velocity_axis_to(\"km/s\")[71]) # 6256.365" ] }, { "cell_type": "markdown", "id": "762c40c0-8219-47b0-adc9-057516bdec08", "metadata": {}, "source": [ "The difference between GCRS (geocentric) and ITRS (topocentric) is very small, a mere 0.110 km/s: 6256.475 vs. 6256.365 km/s" ] }, { "cell_type": "code", "execution_count": 21, "id": "4517728b-5d89-4bab-83c0-d64dda0cc9e7", "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "c3ed9a6732df4b9fb9f141b62eb8672e", "version_major": 2, "version_minor": 0 }, "text/plain": [ "VBox(children=(HBox(children=(Button(description='Clear All Regions', style=ButtonStyle(), tooltip='Clear all …" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "ta.plot(xaxis_unit=\"km/s\");" ] }, { "cell_type": "markdown", "id": "6891ca8c-f3d7-4863-88a2-972835bf2558", "metadata": {}, "source": [ "Here the original the spectral axis has changed, albeit by a very small change.\n", "\n", "500 Hz or 0.10 km/s" ] }, { "cell_type": "code", "execution_count": 22, "id": "a53b8262-43c8-4904-b7be-bf958937eab0", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$[499.59287,~499.59232,~499.59177,~\\dots,~499.51233,~499.51179,~499.51124] \\; \\mathrm{Hz}$" ], "text/plain": [ "" ] }, "execution_count": 22, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ta.spectral_axis - sa" ] }, { "cell_type": "code", "execution_count": 23, "id": "fbb76db6-d1f9-4f69-a454-b55a040ad9c3", "metadata": {}, "outputs": [], "source": [ "ta.set_convention('radio')" ] }, { "cell_type": "code", "execution_count": 24, "id": "697271dd-f0b0-4436-ad09-34944fca8f1d", "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "9d7e6d78c8574869905406d84ebf5624", "version_major": 2, "version_minor": 0 }, "text/plain": [ "VBox(children=(HBox(children=(Button(description='Clear All Regions', style=ButtonStyle(), tooltip='Clear all …" ] }, "metadata": {}, "output_type": "display_data" }, { "name": "stdout", "output_type": "stream", "text": [ "6128.470305969091 km / s\n" ] } ], "source": [ "ta.plot(xaxis_unit=\"km/s\")\n", "print(ta.velocity_axis_to(\"km/s\")[71])" ] }, { "cell_type": "markdown", "id": "a566e5e6-a7dc-4e77-abf3-ce845cf7ecfe", "metadata": {}, "source": [ "The `target` and `observer` are essential if you want to convert between different frames and conventions, as they\n", "control where on the sky and on earth the observations were taken resp." ] }, { "cell_type": "code", "execution_count": 25, "id": "791ccdc7-0118-465c-9ca0-e589c6b0f1f6", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 25, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ta.target" ] }, { "cell_type": "code", "execution_count": 26, "id": "ba74466a-17ac-478a-8450-e3d269fd9b56", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 26, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ta.observer" ] }, { "cell_type": "markdown", "id": "d1fe4839-a722-44ec-bd1f-04f1e64db7db", "metadata": {}, "source": [ "## Other useful functions\n", "\n", "### Spectral Axis Conversion\n", "\n", "Convert the spectral axis to any units, frame, and convention with \n", "[Spectrum.velocity_axis_to](https://dysh.readthedocs.io/en/latest/reference/modules/dysh.spectra.html#dysh.spectra.spectrum.Spectrum.velocity_axis_to). \n", "\n", "dysh understands some common synonyms like 'heliocentric' for astropy's 'hcrs'. Note this returns an array, and does not modify the spectrum meta data.\n" ] }, { "cell_type": "code", "execution_count": 27, "id": "e9cdf384-3249-489e-9a0e-ee96c40c88d9", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$[6220.2727,~6220.6021,~6220.9315,~\\dots,~6268.6937,~6269.0231,~6269.3525] \\; \\mathrm{\\frac{pc}{Myr}}$" ], "text/plain": [ "\n", " target: \n", " observer to target (computed from above):\n", " radial_velocity=4386.005213365332 km / s\n", " redshift=0.014738742242065506\n", " doppler_rest=1420405000.0 Hz\n", " doppler_convention=radio)\n", " [6220.27269497, 6220.60208987, 6220.93148476, 6221.26087966,\n", " 6221.59027455, 6221.91966945, 6222.24906434, 6222.57845924,\n", " 6222.90785413, 6223.23724902, 6223.56664392, 6223.89603881,\n", " 6224.22543371, 6224.5548286 , 6224.8842235 , 6225.21361839,\n", " 6225.54301328, 6225.87240818, 6226.20180307, 6226.53119797,\n", " 6226.86059286, 6227.18998776, 6227.51938265, 6227.84877754,\n", " 6228.17817244, 6228.50756733, 6228.83696223, 6229.16635712,\n", " 6229.49575202, 6229.82514691, 6230.15454181, 6230.4839367 ,\n", " 6230.81333159, 6231.14272649, 6231.47212138, 6231.80151628,\n", " 6232.13091117, 6232.46030607, 6232.78970096, 6233.11909585,\n", " 6233.44849075, 6233.77788564, 6234.10728054, 6234.43667543,\n", " 6234.76607033, 6235.09546522, 6235.42486012, 6235.75425501,\n", " 6236.0836499 , 6236.4130448 , 6236.74243969, 6237.07183459,\n", " 6237.40122948, 6237.73062438, 6238.06001927, 6238.38941416,\n", " 6238.71880906, 6239.04820395, 6239.37759885, 6239.70699374,\n", " 6240.03638864, 6240.36578353, 6240.69517842, 6241.02457332,\n", " 6241.35396821, 6241.68336311, 6242.012758 , 6242.3421529 ,\n", " 6242.67154779, 6243.00094269, 6243.33033758, 6243.65973247,\n", " 6243.98912737, 6244.31852226, 6244.64791716, 6244.97731205,\n", " 6245.30670695, 6245.63610184, 6245.96549673, 6246.29489163,\n", " 6246.62428652, 6246.95368142, 6247.28307631, 6247.61247121,\n", " 6247.9418661 , 6248.27126099, 6248.60065589, 6248.93005078,\n", " 6249.25944568, 6249.58884057, 6249.91823547, 6250.24763036,\n", " 6250.57702526, 6250.90642015, 6251.23581504, 6251.56520994,\n", " 6251.89460483, 6252.22399973, 6252.55339462, 6252.88278952,\n", " 6253.21218441, 6253.5415793 , 6253.8709742 , 6254.20036909,\n", " 6254.52976399, 6254.85915888, 6255.18855378, 6255.51794867,\n", " 6255.84734356, 6256.17673846, 6256.50613335, 6256.83552825,\n", " 6257.16492314, 6257.49431804, 6257.82371293, 6258.15310783,\n", " 6258.48250272, 6258.81189761, 6259.14129251, 6259.4706874 ,\n", " 6259.8000823 , 6260.12947719, 6260.45887209, 6260.78826698,\n", " 6261.11766187, 6261.44705677, 6261.77645166, 6262.10584656,\n", " 6262.43524145, 6262.76463635, 6263.09403124, 6263.42342614,\n", " 6263.75282103, 6264.08221592, 6264.41161082, 6264.74100571,\n", " 6265.07040061, 6265.3997955 , 6265.7291904 , 6266.05858529,\n", " 6266.38798018, 6266.71737508, 6267.04676997, 6267.37616487,\n", " 6267.70555976, 6268.03495466, 6268.36434955, 6268.69374444,\n", " 6269.02313934, 6269.35253423] pc / Myr>" ] }, "execution_count": 27, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ta.velocity_axis_to(unit=\"pc/Myr\", toframe='heliocentric', doppler_convention='radio')" ] }, { "cell_type": "markdown", "id": "021b3360-468e-4c4c-8114-896ab106a8d2", "metadata": {}, "source": [ "### Spectral Shift\n", "\n", "Shift a spectrum in place to a given radial velocity or redshift with \n", "[Spectrum.shift_spectrum_to](https://specutils.readthedocs.io/en/stable/api/specutils.Spectrum.html#specutils.Spectrum.shift_spectrum_to)\n" ] }, { "cell_type": "code", "execution_count": 28, "id": "0a6fc043-c962-487d-a17b-9442f0245478", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "6128.470305969091 km / s\n" ] } ], "source": [ "print(ta.velocity_axis_to(\"km/s\")[71])" ] }, { "cell_type": "code", "execution_count": 29, "id": "df9e7b66-daef-47c8-96d1-8837d083cd49", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "before shift 1391476883.6124058 Hz...\n", "after shift 1412098366.9399488 Hz...\n" ] }, { "name": "stderr", "output_type": "stream", "text": [ "WARNING: FITSFixedWarning: 'obsfix' made the change 'Set OBSGEO-L to 154.294870 from OBSGEO-[XYZ].\n", "Set OBSGEO-B to -1.108897 from OBSGEO-[XYZ].\n", "Set OBSGEO-H to -6378059.370 from OBSGEO-[XYZ]'. [astropy.wcs.wcs]\n" ] } ], "source": [ "print(f\"before shift {ta.spectral_axis[0]}...\")\n", "ta.shift_spectrum_to(radial_velocity=0*u.km/u.s)\n", "print(f\"after shift {ta.spectral_axis[0]}...\")" ] }, { "cell_type": "markdown", "id": "01c5ca68-cc88-48e6-878e-70eb62541718", "metadata": {}, "source": [ "### Spectral Axis in Wavelengths\n", "The default is angstrom, use [Quantity.to](https://docs.astropy.org/en/stable/api/astropy.units.Quantity.html#astropy.units.Quantity.to) to convert to other units." ] }, { "cell_type": "code", "execution_count": 30, "id": "6281df43-02b3-49cf-a348-426bee9e4289", "metadata": {}, "outputs": [ { "data": { "text/latex": [ "$[21.230281,~21.230305,~21.230328,~\\dots,~21.233704,~21.233728,~21.233751] \\; \\mathrm{cm}$" ], "text/plain": [ "\n", " target: \n", " observer to target (computed from above):\n", " radial_velocity=-7.723204298315522e-09 km / s\n", " redshift=-2.5757174171303632e-14\n", " doppler_rest=1420405000.0 Hz\n", " doppler_convention=radio)\n", " [21.23028147, 21.23030476, 21.23032804, 21.23035132, 21.2303746 ,\n", " 21.23039788, 21.23042116, 21.23044444, 21.23046772, 21.230491 ,\n", " 21.23051429, 21.23053757, 21.23056085, 21.23058413, 21.23060741,\n", " 21.23063069, 21.23065398, 21.23067726, 21.23070054, 21.23072382,\n", " 21.2307471 , 21.23077039, 21.23079367, 21.23081695, 21.23084023,\n", " 21.23086351, 21.2308868 , 21.23091008, 21.23093336, 21.23095664,\n", " 21.23097992, 21.23100321, 21.23102649, 21.23104977, 21.23107306,\n", " 21.23109634, 21.23111962, 21.2311429 , 21.23116619, 21.23118947,\n", " 21.23121275, 21.23123603, 21.23125932, 21.2312826 , 21.23130588,\n", " 21.23132917, 21.23135245, 21.23137573, 21.23139902, 21.2314223 ,\n", " 21.23144558, 21.23146887, 21.23149215, 21.23151543, 21.23153872,\n", " 21.231562 , 21.23158529, 21.23160857, 21.23163185, 21.23165514,\n", " 21.23167842, 21.23170171, 21.23172499, 21.23174827, 21.23177156,\n", " 21.23179484, 21.23181813, 21.23184141, 21.23186469, 21.23188798,\n", " 21.23191126, 21.23193455, 21.23195783, 21.23198112, 21.2320044 ,\n", " 21.23202769, 21.23205097, 21.23207426, 21.23209754, 21.23212083,\n", " 21.23214411, 21.2321674 , 21.23219068, 21.23221397, 21.23223725,\n", " 21.23226054, 21.23228382, 21.23230711, 21.23233039, 21.23235368,\n", " 21.23237696, 21.23240025, 21.23242353, 21.23244682, 21.23247011,\n", " 21.23249339, 21.23251668, 21.23253996, 21.23256325, 21.23258654,\n", " 21.23260982, 21.23263311, 21.23265639, 21.23267968, 21.23270297,\n", " 21.23272625, 21.23274954, 21.23277282, 21.23279611, 21.2328194 ,\n", " 21.23284268, 21.23286597, 21.23288926, 21.23291254, 21.23293583,\n", " 21.23295912, 21.2329824 , 21.23300569, 21.23302898, 21.23305226,\n", " 21.23307555, 21.23309884, 21.23312213, 21.23314541, 21.2331687 ,\n", " 21.23319199, 21.23321527, 21.23323856, 21.23326185, 21.23328514,\n", " 21.23330842, 21.23333171, 21.233355 , 21.23337829, 21.23340158,\n", " 21.23342486, 21.23344815, 21.23347144, 21.23349473, 21.23351801,\n", " 21.2335413 , 21.23356459, 21.23358788, 21.23361117, 21.23363446,\n", " 21.23365774, 21.23368103, 21.23370432, 21.23372761, 21.2337509 ] cm>" ] }, "execution_count": 30, "metadata": {}, "output_type": "execute_result" } ], "source": [ "ta.wavelength.to('cm')" ] }, { "cell_type": "markdown", "id": "9df34e85-0840-4fdf-babd-ad966a51bb8a", "metadata": {}, "source": [ "### Plot in Channels\n", "\n", "Also, you can plot the x-axis in channel units" ] }, { "cell_type": "code", "execution_count": 31, "id": "aee63b31-050f-4c5a-ae50-a80f599c7d7a", "metadata": {}, "outputs": [ { "data": { "application/vnd.jupyter.widget-view+json": { "model_id": "d601bfe586df4f8e9d2f1730f1bef0c1", "version_major": 2, "version_minor": 0 }, "text/plain": [ "VBox(children=(HBox(children=(Button(description='Clear All Regions', style=ButtonStyle(), tooltip='Clear all …" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "ta.plot(xaxis_unit='chan');" ] }, { "cell_type": "markdown", "id": "2532b19d-4140-4782-b92d-d84bd891e74e", "metadata": {}, "source": [ "## Final Stats\n", "\n", "Finally, at the end we compute some statistics over a spectrum, merely as a checksum if the notebook is reproducible." ] }, { "cell_type": "code", "execution_count": 32, "id": "577deb46-0400-4c9a-a4a8-4dba704e6d20", "metadata": {}, "outputs": [ { "name": "stderr", "output_type": "stream", "text": [ "20:22:13.670 I rms is OK \n" ] } ], "source": [ "ta.check_stats(0.07586939 * u.K)" ] } ], "metadata": { "kernelspec": { "display_name": "Python 3 (ipykernel)", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3.13.9" }, "toc": { "base_numbering": 0 } }, "nbformat": 4, "nbformat_minor": 5 }