diff --git a/astropy_helpers b/astropy_helpers index 79c3c01a..468f0e7a 160000 --- a/astropy_helpers +++ b/astropy_helpers @@ -1 +1 @@ -Subproject commit 79c3c01af3c436963f275e0f9e7dd3ac0625d70d +Subproject commit 468f0e7a1e5f05e180ebfc1f16caa66cd94d5770 diff --git a/examples/Scintillometry-walk-through.ipynb b/examples/Scintillometry-walk-through.ipynb new file mode 100644 index 00000000..5b67642a --- /dev/null +++ b/examples/Scintillometry-walk-through.ipynb @@ -0,0 +1,389 @@ +{ + "cells": [ + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "# Baseband and Scintillometry Walk Through" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Welcome to the Scintillometry package walk through notebook. Here you will learn the knowledge of pulsar and the basic usage of baseband and scintillometry python packages." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Outline\n", + "* Introducation\n", + "* Read baseband data using baseband pacakge\n", + "* Dedispersion\n", + "* Integration\n", + " * Waterfall plot\n", + " * Folding\n", + " * Polarizations \n", + "* Output" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Sample data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "The sample data used in this notebook are stored at:\n" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Introduction" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Baseband data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Radio telescope receives the electromagnetic waves(EM waves) from far away sources. The EM wave field creates voltage in the telescope antenna. These voltage signals are amplified, filtered and than digitized by the backend electronics. Baseband data are fourier transformed voltage signals. Thus they are complex numbers. " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### Structure of baseband data" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### Polarizations \n", + "\n", + "The EM waves has polarizations. The radio antennas are generally designed to receive the polarized signal by two orthogonal linear polarization channels X and Y. From these two polarization channels, the [Stokes parameters]( https://en.wikipedia.org/wiki/Stokes_parameters) can be calculated. " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Pulsars " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Useful links:
\n", + "For pulsar ephemeris (Pulsar motion parameters)
\n", + "http://www.atnf.csiro.au/research/pulsar/psrcat/
\n", + "For pulse profile
\n", + "http://www.epta.eu.org/epndb/" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Pulsar data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Pulsar data are periodical signals " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "### Pulsar data process" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Baseband --> dedispersion --> Resample --> fold --> Dynamic spectrum --> Secondary Spectrum" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "#### Required python package\n", + "1. [Baseband](https://github.com/mhvk/baseband) \n", + "2. [Numpy](http://www.numpy.org/) \n", + "3. [Astropy](http://www.astropy.org/) \n", + "4. [Scintillometry](https://github.com/mhvk/scintillometry) \n", + "#### Optional python package \n", + "5. [PINT](https://github.com/nanograv/PINT) " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Read baseband data" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "1. Decide your data types\n", + " * Supported data types\n", + " * [VDIF](https://baseband.readthedocs.io/en/stable/vdif/index.html)\n", + " * [MARK 5B](https://baseband.readthedocs.io/en/stable/mark5b/index.html)\n", + " * [MARK 4](https://baseband.readthedocs.io/en/stable/mark4/index.html)\n", + " * [DADA](https://baseband.readthedocs.io/en/stable/dada/index.html)\n", + " * [GUPPI](https://baseband.readthedocs.io/en/stable/guppi/index.html)\n", + " * [GSB](https://baseband.readthedocs.io/en/stable/gsb/index.html)\n", + "\n", + "\n", + "Note: There is a detailed example code in the link to each file types." + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "2. Import the necessary module. For instance, we use GUPPI file type" + ] + }, + { + "cell_type": "code", + "execution_count": 2, + "metadata": {}, + "outputs": [], + "source": [ + "import numpy as np \n", + "import astropy.units as u # import astropy units module, this helps you on the physical quantity unit convers\n", + "from baseband import guppi # import the guppi baseband reader \n", + "import os # This package will help you navigating the file system. " + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "3. Create the file handler for your files. \n", + "
3.1 Read a single file." + ] + }, + { + "cell_type": "code", + "execution_count": 5, + "metadata": {}, + "outputs": [], + "source": [ + "# I think the baseband part should go to the baseband package.\n", + "data_dir = '/home/luo/Projects/baseband/baseband/data'\n", + "data_file = os.path.join(data_dir, 'sample_puppi.raw')\n", + "fh = guppi.open(data_file, 'rb') # 'rb' means read binary file, since the data file is in binary format. " + ] + }, + { + "cell_type": "code", + "execution_count": 12, + "metadata": {}, + "outputs": [ + { + "ename": "UnicodeDecodeError", + "evalue": "'ascii' codec can't decode byte 0xf6 in position 0: ordinal not in range(128)", + "output_type": "error", + "traceback": [ + "\u001b[0;31m---------------------------------------------------------------------------\u001b[0m", + "\u001b[0;31mUnicodeDecodeError\u001b[0m Traceback (most recent call last)", + "\u001b[0;32m\u001b[0m in \u001b[0;36m\u001b[0;34m()\u001b[0m\n\u001b[0;32m----> 1\u001b[0;31m \u001b[0mfh\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mread_header\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m", + "\u001b[0;32m~/.local/lib/python3.6/site-packages/baseband-2.0.0-py3.6.egg/baseband/guppi/base.py\u001b[0m in \u001b[0;36mread_header\u001b[0;34m(self)\u001b[0m\n\u001b[1;32m 98\u001b[0m \u001b[0mheader\u001b[0m \u001b[0;34m:\u001b[0m\u001b[0;31m \u001b[0m\u001b[0;31m`\u001b[0m\u001b[0;34m~\u001b[0m\u001b[0mbaseband\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mguppi\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mGUPPIHeader\u001b[0m\u001b[0;31m`\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 99\u001b[0m \"\"\"\n\u001b[0;32m--> 100\u001b[0;31m \u001b[0;32mreturn\u001b[0m \u001b[0mGUPPIHeader\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfromfile\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mfh_raw\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 101\u001b[0m \u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 102\u001b[0m \u001b[0;32mdef\u001b[0m \u001b[0mread_frame\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0mself\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mmemmap\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mTrue\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0mverify\u001b[0m\u001b[0;34m=\u001b[0m\u001b[0;32mTrue\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;32m~/.local/lib/python3.6/site-packages/baseband-2.0.0-py3.6.egg/baseband/guppi/header.py\u001b[0m in \u001b[0;36mfromfile\u001b[0;34m(cls, fh, verify)\u001b[0m\n\u001b[1;32m 116\u001b[0m \u001b[0mline\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0;34m'========='\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 117\u001b[0m \u001b[0;32mwhile\u001b[0m \u001b[0mline\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;36m8\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;32min\u001b[0m \u001b[0;34m(\u001b[0m\u001b[0;34m'='\u001b[0m\u001b[0;34m,\u001b[0m \u001b[0;34m' '\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0;32m--> 118\u001b[0;31m \u001b[0mline\u001b[0m \u001b[0;34m=\u001b[0m \u001b[0mfh\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mread\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;36m80\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m.\u001b[0m\u001b[0mdecode\u001b[0m\u001b[0;34m(\u001b[0m\u001b[0;34m'ascii'\u001b[0m\u001b[0;34m)\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[0m\u001b[1;32m 119\u001b[0m \u001b[0;32mif\u001b[0m \u001b[0mline\u001b[0m\u001b[0;34m[\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;36m3\u001b[0m\u001b[0;34m]\u001b[0m \u001b[0;34m==\u001b[0m \u001b[0;34m'END'\u001b[0m\u001b[0;34m:\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n\u001b[1;32m 120\u001b[0m \u001b[0;32mbreak\u001b[0m\u001b[0;34m\u001b[0m\u001b[0m\n", + "\u001b[0;31mUnicodeDecodeError\u001b[0m: 'ascii' codec can't decode byte 0xf6 in position 0: ordinal not in range(128)" + ] + } + ], + "source": [ + "fh.read_header()" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "
3.2 Read mulitple files together. \n", + "Note: Those files should be contiguous and sorted by time. " + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "4. Read and exam the header information" + ] + }, + { + "cell_type": "code", + "execution_count": 11, + "metadata": {}, + "outputs": [ + { + "data": { + "text/plain": [ + "GUPPIFileReader(fh_raw=<_io.BufferedReader name='/home/luo/Projects/baseband/baseband/data/sample_puppi.raw'>)" + ] + }, + "execution_count": 11, + "metadata": {}, + "output_type": "execute_result" + } + ], + "source": [ + "fh." + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [ + "\"\"\"Test of the scintillometry dedispersion routines on PSR B1937+21\"\"\"\n", + "import numpy as np\n", + "import astropy.units as u\n", + "from baseband import guppi\n", + "from scintillometry.dispersion import Dedisperse\n", + "from scintillometry.functions import Square\n", + "from scintillometry.fourier import get_fft_maker\n", + "import matplotlib.pylab as plt\n", + "\n", + "\n", + "FFT = get_fft_maker('pyfftw', flags=['FFTW_ESTIMATE'], threads=4)\n", + "\n", + "dm = 71.02227638 * u.pc / u.cm**3\n", + "f0 = (641.928224582360 + 2.78557046583437451/60.) * u.Hz # from polyco\n", + "\n", + "fh = guppi.open('puppi_58245_B1937+21_0799.0010.raw')\n", + "frequency = (fh.header0['OBSFREQ']\n", + " - fh.header0['OBSBW'] / 2. # Down to bottom end\n", + " + fh.header0['CHAN_BW'] / 2. # Up to center of channel 0\n", + " + fh.header0['CHAN_BW'] * np.arange(4)) * u.MHz\n", + "\n", + "dedisperse = Dedisperse(fh, dm, frequency=frequency, sideband=1,\n", + " FFT=FFT)\n", + "\n", + "square = Square(dedisperse)\n", + "\n", + "nbin = 64\n", + "lc = np.zeros(nbin)\n", + "for i in range(10):\n", + " print('At ', square.tell(unit=u.s))\n", + " pre = square.tell()\n", + " p = square.read(3125000).sum((1, 2))\n", + " post = square.tell()\n", + " ph = (np.arange(pre, post) / fh.sample_rate * f0).to_value(u.one)\n", + " ph = (ph - 0.35) % 1.\n", + " lc += np.bincount((ph * nbin).astype(int), weights=p)\n", + "\n", + "plt.ion()\n", + "plt.plot(np.arange(nbin)/nbin, lc)" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "Useful links:
\n", + "    Baseband documentation: https://baseband.readthedocs.io/en/stable/ \n", + "
\n", + "    Baseband github: https://github.com/mhvk/baseband" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Dedispersion" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Integration" + ] + }, + { + "cell_type": "markdown", + "metadata": {}, + "source": [ + "## Output" + ] + }, + { + "cell_type": "code", + "execution_count": null, + "metadata": {}, + "outputs": [], + "source": [] + } + ], + "metadata": { + "kernelspec": { + "display_name": "Python 3", + "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.6.5" + } + }, + "nbformat": 4, + "nbformat_minor": 2 +}